Nothing
## =========================================================================================
##
## PUBLIC S3 METHODS FOR PLNmixturefit
##
## =========================================================================================
## Auxiliary functions to check the given class of an objet
isPLNmixturefit <- function(Robject) {inherits(Robject, "PLNmixturefit")}
#' Mixture visualization of a [`PLNmixturefit`] object
#'
#' Represent the result of the clustering either by coloring the individual in a two-dimension PCA factor map,
#' or by representing the expected matrix of count reorder according to the clustering.
#'
#' @name plot.PLNmixturefit
#'
#' @param x an R6 object with class [`PLNmixturefit`]
#' @param type character for the type of plot, either "pca", for or "matrix". Default is `"pca"`.
#' @param main character. A title for the plot. If NULL (the default), an hopefully appropriate title will be used.
#' @param plot logical. Should the plot be displayed or sent back as [`ggplot`] object
#' @param ... Not used (S3 compatibility).
#'
#' @return a [`ggplot`] graphic
#' @examples
#' data(trichoptera)
#' trichoptera <- prepare_data(trichoptera$Abundance, trichoptera$Covariate)
#' myPLN <- PLNmixture(Abundance ~ 1 + offset(log(Offset)),
#' data = trichoptera, control = PLNmixture_param(smoothing = "none")) %>% getBestModel()
#' \dontrun{
#' plot(myPLN, "pca")
#' plot(myPLN, "matrix")
#' }
#' @export
plot.PLNmixturefit <-
function(x,
type = c("pca", "matrix"),
main = NULL,
plot = TRUE, ...) {
if (is.null(main))
main <- switch(match.arg(type),
"pca" = "Clustering labels in Individual Factor Map",
"matrix" = "Expected counts reorder by clustering")
p <- switch(match.arg(type),
"pca" = x$plot_clustering_pca(main = main, plot = FALSE),
"matrix" = x$plot_clustering_data(main = main, plot = FALSE))
if (plot) print(p)
invisible(p)
}
#' Prediction for a [`PLNmixturefit`] object
#'
#' Predict either posterior probabilities for each group or latent positions based on new samples
#'
#' @param object an R6 object with class [`PLNmixturefit`]
#' @param newdata A data frame in which to look for variables, offsets and counts with which to predict.
#' @param type The type of prediction required. The default `posterior` are posterior probabilities for each group ,
#' `response` is the group with maximal posterior probability and `latent` is the averaged latent in the latent space,
#' with weights equal to the posterior probabilities.
#' @param prior User-specified prior group probabilities in the new data. The default uses a uniform prior.
#' @param control a list-like structure for controlling the fit. See [PLNmixture_param()] for details.
#' @param ... additional parameters for S3 compatibility. Not used
#' @return A matrix of posterior probabilities for each group (if type = "posterior"), a matrix of (average) position in the
#' latent space (if type = "position") or a vector of predicted groups (if type = "response").
#' @export
#' @examples
#' data(trichoptera)
#' trichoptera <- prepare_data(trichoptera$Abundance, trichoptera$Covariate)
#' myPLN <- PLNmixture(Abundance ~ 1 + offset(log(Offset)),
#' data = trichoptera, control = PLNmixture_param(smoothing = "none")) %>% getBestModel()
#' predict(myPLN, trichoptera, "posterior")
#' predict(myPLN, trichoptera, "position")
#' predict(myPLN, trichoptera, "response")
predict.PLNmixturefit <-
function(object, newdata,
type = c("posterior", "response", "position"),
prior = matrix(rep(1/object$k, object$k), nrow(newdata), object$k, byrow = TRUE),
control = PLNmixture_param(), ...) {
stopifnot(isPLNmixturefit(object))
object$predict(newdata, type, prior, control, parent.frame())
}
#' Extract model coefficients
#'
#' @description Extracts model coefficients from objects returned by [PLN()] and its variants
#'
#' @name coef.PLNmixturefit
#'
#' @param object an R6 object with class [`PLNmixturefit`]
#' @param type type of parameter that should be extracted. Either "main" (default) for \deqn{\Theta},
#' "means" for \deqn{\mu}, "mixture" for \deqn{\pi} or "covariance" for \deqn{\Sigma}
#' @param ... additional parameters for S3 compatibility. Not used
#' @return A matrix of coefficients extracted from the PLNfit model.
#'
#' @seealso [sigma.PLNmixturefit()]
#'
#' @export
#' @examples
#' data(trichoptera)
#' trichoptera <- prepare_data(trichoptera$Abundance, trichoptera$Covariate)
#' myPLN <- PLNmixture(Abundance ~ 1 + offset(log(Offset)),
#' data = trichoptera, control = PLNmixture_param(smoothing = "none")) %>% getBestModel()
#' coef(myPLN) ## Theta - empty here
#' coef(myPLN, type = "mixture") ## pi
#' coef(myPLN, type = "means") ## mu
#' coef(myPLN, type = "covariance") ## Sigma
coef.PLNmixturefit <- function(object, type = c("main", "means", "covariance", "mixture"), ...) {
stopifnot(isPLNmixturefit(object))
switch(match.arg(type),
main = object$model_par$Theta,
means = object$group_means,
mixture = object$mixtureParam,
covariance = object$model_par$Sigma)
}
#' Extracts model fitted values from objects returned by [PLNmixture()] and its variants
#'
#' @name fitted.PLNmixturefit
#'
#' @inheritParams coef.PLNmixturefit
#' @return A matrix of Fitted values extracted from the object object.
#'
#' @export
fitted.PLNmixturefit <- function(object, ...) {
stopifnot(isPLNmixturefit(object))
object$fitted
}
#' Extract variance-covariance of residuals 'Sigma'
#'
#' @name sigma.PLNmixturefit
#' @description Extract the variance-covariance matrix of the residuals, usually noted \deqn{\Sigma} in PLN models. This captures the correlation between the species in the latent space.
#' or PLNmixture, it is a weighted mean of the variance-covariance matrices of each component.
#'
#' @inheritParams coef.PLNmixturefit
#'
#' @return A semi definite positive matrix of size p, assuming there are p species in the model.
#'
#' @export
#'
#' @seealso [coef.PLNmixturefit()] for other ways to access \deqn{\Sigma}.
#'
#' @importFrom stats sigma
#' @examples
#' data(trichoptera)
#' trichoptera <- prepare_data(trichoptera$Abundance, trichoptera$Covariate)
#' myPLN <- PLNmixture(Abundance ~ 1 + offset(log(Offset)),
#' data = trichoptera, control = PLNmixture_param(smoothing = "none")) %>% getBestModel()
#' sigma(myPLN) ## Sigma
sigma.PLNmixturefit <- function(object, ...) {
stopifnot(isPLNmixturefit(object))
object$model_par$Sigma
}
#' @describeIn standard_error Component-wise standard errors of B in [`PLNmixturefit`] (not implemented yet)
#' @export
standard_error.PLNmixturefit <- function(object, type = c("variational", "jackknife", "sandwich"), parameter = c("B", "Omega")) {
par <- match.arg(parameter)
if (par == "Omega")
stop("Omega is not estimated as such in PLNmixture models")
if (par == "B") {
warning("Standard error of B is not implemented yet for PLNmixture models")
theta_sd <- coef.PLNmixturefit(object)
theta_sd[ , ] <- NA
theta_sd
}
}
Any scripts or data that you put into this service are public.
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.