Description Usage Arguments Value About type About spd About solver About detap Examples
coxmeg
returns estimates of the variance component, the HRs and p-values for the predictors.
1 2 3 4 |
outcome |
A matrix contains time (first column) and status (second column). The status is a binary variable (1 for events / 0 for censored). |
corr |
A relatedness matrix. Can be a matrix or a 'dgCMatrix' class in the Matrix package. Must be symmetric positive definite or symmetric positive semidefinite. |
type |
A string indicating the sparsity structure of the relatedness matrix. Should be 'bd' (block diagonal), 'sparse', or 'dense'. See details. |
X |
An optional matrix of the preidctors with fixed effects. Can be quantitative or binary values. Categorical variables need to be converted to dummy variables. Each row is a sample, and the predictors are columns. |
FID |
An optional string vector of family ID. If provided, the data will be reordered according to the family ID. |
eps |
An optional positive scalar indicating the tolerance in the optimization algorithm. Default is 1e-6. |
min_tau |
An optional positive scalar indicating the lower bound in the optimization algorithm for the variance component |
max_tau |
An optional positive scalar indicating the upper bound in the optimization algorithm for the variance component |
order |
An optional integer starting from 0. Only valid when |
detap |
An optional string indicating whether to use approximation for log-determinant. Can be 'exact', 'diagonal' or 'slq'. Default is NULL, which lets the function select a method based on 'type' and other information. See details. |
opt |
An optional logical scalar for the Optimization algorithm for tau. Can have the following values: 'bobyqa', 'Brent' or 'NM'. Default is 'bobyqa'. |
solver |
An optional bianry value that can be either 1 (Cholesky Decomposition using RcppEigen), 2 (PCG) or 3 (Cholesky Decomposition using Matrix). Default is NULL, which lets the function select a solver. See details. |
spd |
An optional logical value indicating whether the relatedness matrix is symmetric positive definite. Default is TRUE. See details. |
verbose |
An optional logical scalar indicating whether to print additional messages. Default is TRUE. |
mc |
An optional integer scalar specifying the number of Monte Carlo samples used for approximating the log-determinant. Only valid when |
beta: The estimated coefficient for each predictor in X.
HR: The estimated HR for each predictor in X.
sd_beta: The estimated standard error of beta.
p: The p-value.
iter: The number of iterations until convergence.
tau: The estimated variance component.
int_ll: The marginal likelihood (-2*log(lik)) of tau evaluated at the estimate of tau.
rank: The rank of the relatedness matrix.
nsam: Actual sample size.
type
'bd' is used for a block-diagonal relatedness matrix, or a sparse matrix the inverse of which is also sparse. 'sparse' is used for a general sparse relatedness matrix the inverse of which is not sparse.
spd
When spd=TRUE
, the relatedness matrix is treated as SPD. If the matrix is SPSD or not sure, use spd=FALSE
.
solver
When solver=1,3
/solver=2
, Cholesky decompositon/PCG is used to solve the linear system. When solver=3
, the solve function in the Matrix package is used, and when solver=1
, it uses RcppEigen:LDLT to solve linear systems. When type='dense'
, it is recommended to set solver=2
to have better computational performance.
detap
When detap='exact'
, the exact log-determinant is computed for estimating the variance component. Specifying detap='diagonal'
uses diagonal approximation, and is only effective for a sparse relatedness matrix. Specifying detap='slq'
uses stochastic lanczos quadrature approximation.
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 | library(Matrix)
library(MASS)
library(coxmeg)
## simulate a block-diagonal relatedness matrix
tau_var <- 0.2
n_f <- 100
mat_list <- list()
size <- rep(10,n_f)
offd <- 0.5
for(i in 1:n_f)
{
mat_list[[i]] <- matrix(offd,size[i],size[i])
diag(mat_list[[i]]) <- 1
}
sigma <- as.matrix(bdiag(mat_list))
n <- nrow(sigma)
## simulate random effects and outcomes
x <- mvrnorm(1, rep(0,n), tau_var*sigma)
myrates <- exp(x-1)
y <- rexp(n, rate = myrates)
cen <- rexp(n, rate = 0.02 )
ycen <- pmin(y, cen)
outcome <- cbind(ycen,as.numeric(y <= cen))
## fit a Cox mixed-effects model
re = coxmeg(outcome,sigma,type='bd',detap='diagonal',order=1)
re
|
Loading required package: Rcpp
Remove 0 subjects censored before the first failure.
There is/are 0 predictors. The sample size included is 1000.
The relatedness matrix is treated as sparse.
The relatedness matrix is inverted.
The method for computing the determinant is 'diagonal'.
Solver: Cholesky decomposition (RcppEigen=TRUE).
$beta
NULL
$HR
NULL
$sd_beta
NULL
$p
NULL
$tau
[1] 0.08348784
$iter
[1] 21
$rank
[1] 1000
$nsam
[1] 1000
$int_ll
[1] 11582.22
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.