cisAssoc {GGtools} | R Documentation |
test for variant-expression associations in cis, using VCF and SummarizedExperiment representations
cisAssoc(summex, vcf.tf, rhs = ~1, nperm = 3, cisradius = 50000, genome = "hg19", assayind = 1, lbmaf = 1e-06, dropUnivHet = TRUE, doEsts=FALSE) data(lgeu) # obtains an example SummarizedExperiment
summex |
instance of |
vcf.tf |
instance of |
rhs |
formula ‘right hand side’ for adjustments to be made
as |
nperm |
number of permutations to be used for plug-in FDR computation |
cisradius |
distance in bp around each gene body to be searched for SNP association |
genome |
tag suitable for use in GenomeInfoDb structures |
assayind |
index of |
lbmaf |
lower bound on MAF of SNP to use |
dropUnivHet |
logical, if TRUE, will check for columns of SnpMatrix instance that possess no values other than "NA" and "A/B". See http://www.biostars.org/p/117155/#117270 |
doEsts |
logical, if TRUE, will run |
snp.rhs.tests
is the
workhorse for statistical modeling. VCF content is
transformed to the byte-code (which allows for uncertain imputation)
and used in fast testing.
a GRanges-class
instance
with mcols including chisq, permScore...
seqlevelsStyle for summex and vcf.tf content must agree
VJ Carey <stvjd@channing.harvard.edu>
data(lgeu) # small excerpt from GEUVADIS FPKM tf20 = TabixFile(system.file("vcf/c20exch.vcf.gz", package="GGtools")) if (require(VariantAnnotation)) scanVcfHeader(tf20) lgeue = clipPCs(lgeu[,which(lgeu$popcode=="CEU")], 1:2) set.seed(4321) litc = cisAssoc(lgeue, tf20, nperm=2, lbmaf=.05, cisradius=50000) litc2 = cisAssoc(lgeue, tf20, nperm=2, lbmaf=.05, cisradius=50000, doEsts=TRUE) summary(litc$chisq) litc$pifdr = pifdr(litc$chisq, c(litc$permScore_1, litc$permScore_2)) litc[which(litc$pifdr < .01)]