kmFit_eQTL: eQTL linear mixed effects models with kinship

In BIGslu/kimma: Kinship In Mixed Model Analysis of RNA-seq

eQTL linear mixed effects models with kinship

Description

Run lmerel and corresponding lm or lme without kinship of gene expression in RNA-seq data against associated genotypes

Usage

kmFit_eQTL(
  dat_snp = NULL,
  dat_map = NULL,
  geneID = "gene",
  genotypeID = "snpID",
  dat = NULL,
  kin = NULL,
  patientID = "ptID",
  libraryID = "libID",
  counts = NULL,
  meta = NULL,
  genes = NULL,
  weights = NULL,
  subset_var = NULL,
  subset_lvl = NULL,
  subset_genes = NULL,
  model,
  use_weights = FALSE,
  run_lm = FALSE,
  run_lme = FALSE,
  run_lmerel = FALSE,
  metrics = FALSE,
  run_contrast = FALSE,
  contrast_var = NULL,
  processors = NULL,
  p_method = "BH",
  dat.snp = NULL,
  dat.map = NULL,
  run.lmekin = NULL,
  subset.var = NULL,
  subset.lvl = NULL,
  subset.genes = NULL,
  use.weights = FALSE,
  run.lm = FALSE,
  run.lme = FALSE,
  run.lmerel = FALSE,
  run.contrast = FALSE,
  contrast.var = NULL,
  p.method = NULL
)

Arguments

dat_snp	Data frame containing 0,1,2 genotypes. Can be numeric or character/factor depending on desired analysis. If running contrasdts, must be character/factor for pairwise comparisons to run. Rows are donors and columns are SNP.
dat_map	Data frame mapping geneID to genotypeID. If provided, genotypes will be modeled against only the genes they are paired with in the map. If not provided, genotypes will be modeled against all genes in dat$E / counts
geneID	Character variable name to match dat_snp dat with gene identifiers. Values should match rownames in express data (dat$E or counts) Default it "gene"
genotypeID	Character variable name to match dat_snp dat with SNP identifiers. Default it "snpID"
dat	EList object output by voom( ). Contains counts (dat$E), meta (dat$targets), and genes (dat$genes).
kin	Matrix with pairwise kinship values between individuals. Must be numeric with rownames.
patientID	Character of variable name to match dat$targets to kinship row and column names.
libraryID	Character of variable name to match dat$targets to dat$E colnames
counts	Matrix of normalized expression. Rows are genes, columns are libraries.
meta	Matrix or data frame of sample and individual metadata.
genes	Matrix or data frame of gene metadata.
weights	Matrix of data frame of gene specific weights. Usually calculated with limma::voomWithQualityWeights().
subset_var	Character list of variable name(s) to filter data by.
subset_lvl	Character list of variable value(s) or level(s) to filter data to. Must match order of subset.var
subset_genes	Character vector of genes to include in models.
model	Character vector of model starting with ~ Should include (1|patientID) if mixed effects will be run. NOTE THAT FOR EQTL, the model MUST include genotype as a place holder for the genotype data. For example "~ genotype + (1|ptID)" Same as kmFit( )
use_weights	Logical if gene specific weights should be used in model. Default is FALSE
run_lm	Logical if should run lm model without kinship
run_lme	Logical if should run lme model without kinship
run_lmerel	Logical if should run lmerel model with kinship
metrics	Logical if should calculate model fit metrics such as AIC, BIC, R-squared. Default is FALSE
run_contrast	Logical if should run pairwise contrasts. If no matrix provided, all possible pairwise comparisons are completed.
contrast_var	Character vector of variable in model to run contrasts of. Interaction terms must be specified as "var1:var2". If NULL (default), all contrasts for all variables in the model are run
processors	Numeric processors to run in parallel. Default is 2 less than the total available
p_method	Character of FDR adjustment method. Values as in p.adjust( )
dat.snp	Deprecated form of dat_snp
dat.map	Deprecated form of dat.map
run.lmekin	Deprecated. Please use run_lmerel
subset.var	Deprecated form of subset_var
subset.lvl	Deprecated form of subset_lvl
subset.genes	Deprecated form of subset_genes
use.weights	Deprecated form of use_weights
run.lm	Deprecated form of run_lm
run.lme	Deprecated form of run_lme
run.lmerel	Deprecated form of run_lmerel
run.contrast	Deprecated form of run_contrast
contrast.var	Deprecated form of contrast_var
p.method	Deprecated form of p_method

Value

List of data frames including - lm/lme/lmerel: model estimates and significance - *.contrast: model estimates and significance for pairwise contrasts with variables in the original model - *.fit: model fit metrics such as sigma, AIC, BIC, R-squared (optional with metrics paramater) - *.error: error messages for genes that failed model fitting

Examples

#Determine data size
nSNP <- 4
nSamp <- length(unique(example.voom$targets$ptID))
#Create random genotype data
snp <- data.frame(matrix(sample(c(0,1,2), nSNP*nSamp, replace=TRUE),
                         nrow=nSamp)) %>%
  #Add SNP names
  dplyr::mutate(ptID = unique(example.voom$targets$ptID)) %>%
  dplyr::select(ptID, dplyr::everything())
colnames(snp) <- c("ptID", paste0("rs",1:nSNP))

#Create gene to SNP map
 map <- data.frame(snpID = paste0("rs",1:nSNP),
                 gene = rep(example.voom$genes$geneName[1:2], nSNP/2))

#Run SNP interaction model
eQTL_result <- kmFit_eQTL(dat_snp = snp, dat_map = map,
                          dat = example.voom, kin = example.kin,
                          run_lmerel = TRUE, metrics=TRUE,
                          run_contrast = TRUE, contrast_var="virus:genotype",
                          model = "~ virus*genotype + (1|ptID)")

#Run SNP as factors with contrasts
snp_fct <- snp %>%
 dplyr::mutate_if(is.numeric, ~as.factor(.))

eQTL_result2 <- kmFit_eQTL(dat_snp = snp_fct, dat_map = map,
                           dat = example.voom, kin = example.kin,
                           run_lmerel = TRUE, run_contrast = TRUE, metrics=TRUE,
                           model = "~ genotype + (1|ptID)")

BIGslu/kimma documentation built on May 30, 2024, 10:09 p.m.