CJAMP (version 0.1.1)

minusloglik: Minus log-likelihood of copula models.

Description

Function to compute the minus log-likelihood of the joint distribution of Y1 and Y2 given the predictors predictors_Y1 of Y1 and predictors_Y2 of Y2 in a copula model for given parameter values parameters. Implemented are the Clayton and 2-parameter copula. The function assumes quantitative predictors and uses an additive model, i.e. for categorical predictors with more than 2 levels, dummy variables have to be created beforehand. Accordingly, if single nucleotide variants (SNVs) are included as predictors, the computation is based on an additive genetic model if SNVs are provided as 0-1-2 genotypes and on a dominant model if SNVs are provided as 0-1 genotypes.

Usage

minusloglik(copula = "Clayton", Y1 = NULL, Y2 = NULL,
  predictors_Y1 = NULL, predictors_Y2 = NULL, parameters = NULL)

Arguments

copula

String indicating whether the likelihood should be computed under the Clayton ("Clayton") or 2-parameter copula ("2param") model.

Y1

Numeric vector containing the first phenotype.

Y2

Numeric vector containing the second phenotype.

predictors_Y1

Dataframe containing the predictors of Y1 in columns.

predictors_Y2

Dataframe containing predictors of Y2 in columns.

parameters

Named integer vector containing the parameter estimates of the marginal and copula parameters, for which the log-likelihood will be computed. For details and the necessary format of the vector, see the details below.

Value

Minus log-likelihood value (integer).

Details

The number of predictors is not fixed and also different predictors can be supplied for Y1 and Y2. However, for the functioning of the log-likelihood function, the parameters vector has to be supplied in a pre-specified form and formatting:

The vector has to include the copula parameter(s), marginal parameters for Y1, and marginal parameters for Y2 in this order. For example, for the 2-parameter copula with parameters \(\phi, \theta\) and marginal models

$$Y_1 = \alpha_0 + \alpha_1 \cdot X_1 + \alpha_2 \cdot X_2 + \epsilon_1, \ \epsilon_1 \sim N(0,\sigma_1^2),$$ $$Y_1 = \beta_0 + \beta_1 \cdot X_1 + \epsilon_2, \ \epsilon_2 \sim N(0,\sigma_2^2),$$

the parameter vector has be $$(log(\phi), log(\theta-1), log(\sigma_1), log(\sigma_2), \alpha_0, \alpha_2, \alpha_1, \beta_0, \beta_1)^T.$$ \(log(\phi)\) and \(log(\theta-1)\) have to be provided instead of \(\phi, \theta\) and \(log(\sigma_1)\), \(log(\sigma_2)\) instead of \(\sigma_1\), \(\sigma_2\) for computational reasons when the log-likelihood function is optimized in cjamp. As further necessary format, the vector has to be named and the names of the copula parameter(s) has/have to be log_phi and/or log_theta_minus1.

Examples

# NOT RUN {
# Generate genetic data:
set.seed(10)
genodata <- generate_genodata(n_SNV = 20, n_ind = 1000)

# Generate phenotype data:
phenodata <- generate_phenodata_2_copula(genodata = genodata, MAF_cutoff = 1,
                                         prop_causal = 0.5, tau = 0.2,
                                         b1 = 0.3, b2 = 0.3)

# Example 1: Log-likelihood of null model without covariates & genetic effects
estimates <- get_estimates_naive(Y1 = phenodata$Y1, Y2 = phenodata$Y2,
                                 predictors_Y1 = NULL, predictors_Y2 = NULL,
                                 copula_param = "both")
minusloglik(Y1 = phenodata$Y1, Y2 = phenodata$Y2, predictors_Y1 = NULL,
            predictors_Y2 = NULL, parameters = estimates, copula = "2param")

# Example 2: Log-likelihood of null model with covariates, without genetic effects
predictors <- data.frame(X1 = phenodata$X1, X2 = phenodata$X2)
estimates <- get_estimates_naive(Y1 = phenodata$Y1, Y2 = phenodata$Y2,
                                 predictors_Y1 = predictors,
                                 predictors_Y2 = predictors,
                                 copula_param = "both")
minusloglik(Y1 = phenodata$Y1, Y2 = phenodata$Y2, predictors_Y1 = predictors,
            predictors_Y2 = predictors, parameters = estimates, copula = "2param")

# Example 3: Log-likelihood of model with covariates & genetic effect on Y1 only
predictors_Y1 <- data.frame(X1 = phenodata$X1, X2 = phenodata$X2,
                            SNV = genodata$SNV1)
predictors_Y2 <- data.frame(X1 = phenodata$X1, X2 = phenodata$X2)
estimates <- get_estimates_naive(Y1 = phenodata$Y1, Y2 = phenodata$Y2,
                                 predictors_Y1 = predictors_Y1,
                                 predictors_Y2 = predictors_Y2,
                                 copula_param = "both")
minusloglik(Y1 = phenodata$Y1, Y2 = phenodata$Y2, predictors_Y1 = predictors_Y1,
            predictors_Y2 = predictors_Y2, parameters = estimates,
            copula = "2param")

# Example 4: Log-likelihood of model with covariates & genetic effect on Y2 only
predictors_Y1 <- data.frame(X1 = phenodata$X1, X2 = phenodata$X2)
predictors_Y2 <- data.frame(X1 = phenodata$X1, X2 = phenodata$X2,
                            SNV = genodata$SNV1)
estimates <- get_estimates_naive(Y1 = phenodata$Y1, Y2 = phenodata$Y2,
                                 predictors_Y1 = predictors_Y1,
                                 predictors_Y2 = predictors_Y2,
                                 copula_param = "both")
minusloglik(Y1 = phenodata$Y1, Y2 = phenodata$Y2, predictors_Y1 = predictors_Y1,
            predictors_Y2 = predictors_Y2, parameters = estimates,
            copula = "2param")

# Example 5: Log-likelihood of model without covariates, with genetic effects
predictors <- data.frame(SNV = genodata$SNV1)
estimates <- get_estimates_naive(Y1 = phenodata$Y1, Y2 = phenodata$Y2,
                                 predictors_Y1 = predictors,
                                 predictors_Y2 = predictors,
                                 copula_param = "both")
minusloglik(Y1 = phenodata$Y1, Y2 = phenodata$Y2, predictors_Y1 = predictors,
            predictors_Y2 = predictors, parameters = estimates, copula = "2param")

# Example 6: Log-likelihood of model with covariates & genetic effects
predictors <- data.frame(X1 = phenodata$X1, X2 = phenodata$X2, SNV = genodata$SNV1)
estimates <- get_estimates_naive(Y1 = phenodata$Y1, Y2 = phenodata$Y2,
                                 predictors_Y1 = predictors,
                                 predictors_Y2 = predictors,
                                 copula_param = "both")
minusloglik(Y1 = phenodata$Y1, Y2 = phenodata$Y2, predictors_Y1 = predictors,
            predictors_Y2 = predictors, parameters = estimates, copula = "2param")

# Example 7: Log-likelihood of model with covariates & multiple genetic effects
predictors <- data.frame(X1 = phenodata$X1, X2 = phenodata$X2, genodata[, 1:5])
estimates <- get_estimates_naive(Y1 = phenodata$Y1, Y2 = phenodata$Y2,
                                 predictors_Y1 = predictors,
                                 predictors_Y2 = predictors,
                                 copula_param = "both")
minusloglik(Y1 = phenodata$Y1, Y2 = phenodata$Y2, predictors_Y1 = predictors,
            predictors_Y2 = predictors, parameters = estimates, copula = "2param")

# }