coxmeg_plink {coxmeg} | R Documentation |
coxmeg_plink
first estimates the variance component under a null model with only cov if tau is not given, and then analyzing each SNP in the plink files.
coxmeg_plink(pheno, corr, type, bed = NULL, tmp_dir = NULL, cov_file = NULL, tau = NULL, maf = 0.05, min_tau = 1e-04, max_tau = 5, eps = 1e-06, order = NULL, detap = NULL, opt = "bobyqa", score = FALSE, threshold = 0, solver = NULL, spd = TRUE, mc = 100, verbose = TRUE)
pheno |
A string value indicating the file name or the full path of a pheno file. The files must be in the working directory if the full path is not given. The file is in plink pheno format, containing the following four columns, family ID, individual ID, time and status. 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. |
bed |
A optional string value indicating the file name or the full path of a plink bed file (without .bed). The files must be in the working directory if the full path is not given. If not provided, only the variance component will be returned. |
tmp_dir |
A optional directory to store temporary .gds files. The directory needs to be specified when |
cov_file |
An optional string value indicating the file name or the full path of a covariate file. The files must be in the working directory if the full path is not given. Same as the cov file in plink, the first two columns are family ID and individual ID. The covariates are included in the null model for estimating the variance component. The covariates can be quantitative or binary values. Categorical variables need to be converted to dummy variables. |
tau |
An optional positive value for the variance component. If tau is given, the function will skip estimating the variance component, and use the given tau to analyze the SNPs. |
maf |
An optional positive value. All SNPs with MAF<maf in the bed file will not be analyzed. Default is 0.05. |
min_tau |
An optional positive value indicating the lower bound in the optimization algorithm for the variance component tau. Default is 1e-4. |
max_tau |
An optional positive value indicating the upper bound in the optimization algorithm for the variance component tau. Default is 5. |
eps |
An optional positive value indicating the tolerance in the optimization algorithm. Default is 1e-6. |
order |
An optional integer value starting from 0. Only effective 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 string value for the Optimization algorithm for tau. Can have the following values: 'bobyqa', 'Brent' or 'NM'. Default is 'bobyqa'. |
score |
An optional logical value indicating whether to perform a score test. Default is FALSE. |
threshold |
An optional non-negative value. If threshold>0, coxmeg_m will reestimate HRs for those SNPs with a p-value<threshold by first estimating a variant-specific variance component. Default is 0. |
solver |
An optional binary value taking either 1 or 2. Default is 1. See details. |
spd |
An optional logical value indicating whether the relatedness matrix is symmetric positive definite. Default is TRUE. See details. |
mc |
An optional integer value specifying the number of Monte Carlo samples used for approximating the log-determinant. Only valid when |
verbose |
An optional logical value indicating whether to print additional messages. Default is TRUE. |
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 of each SNP.
tau: The estimated variance component.
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.
corr
The subjects in corr
must be in the same order as in the plink fam file.
pheno
-9 for missing values, cov_file
NA for missing values.
The function will create a temporary gds file with approximately the same size as the bed file. The temporary file will be removed when the analysis is done.
spd
When spd=TRUE
, the relatedness matrix is treated as SPD. If the matrix is SPSD or not sure, set 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.
library(Matrix) library(MASS) library(coxmeg) ## build a block-diagonal relatedness matrix n_f <- 600 mat_list <- list() size <- rep(5,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)) ## Estimate variance component under a null model pheno = system.file("extdata", "ex_pheno.txt", package = "coxmeg") cov = system.file("extdata", "ex_cov.txt", package = "coxmeg") bed = system.file("extdata", "example_null.bed", package = "coxmeg") bed = substr(bed,1,nchar(bed)-4) re = coxmeg_plink(pheno,sigma,type='bd',bed=bed,tmp_dir=tempdir(),cov_file=cov, detap='diagonal',order=1) re