rbm.haplos.power {Ravages}

Power of RVAT based on simulations with haplotypes

Description

Computes the power of the tests CAST, WSS and SKAT based on simulations with haplotypes

Usage

  rbm.haplos.power(haplos, freqs, weights = "SKAT", 
                   max.maf.causal = 0.01, maf.filter = max.maf.causal, p.causal = 0.5, 
                   p.protect = 0, h2 = c(0.01, 0.01), prev = c(1, 0.01),
                   normal.approx = TRUE, size = c(500, 500), verbose = TRUE,
                   alpha = 2.5e-6, RVAT = c("CAST", "WSS", "SKAT"),
                   SKAT.method = c("permutations", "theoretical"),
                   simus.haplos = c("freqs", "liability"), 
                   replicates = 1000, rep.by.causal = 50, cores = 10)

Arguments

haplos

A matrix of haplotypes with one row per haplotype and one column per variant

freqs

A matrix of haplotypes frequencies in each group of individuals, only needed if simus.haplos = "freqs"

weights

How to weight rare variants (if "constant", all variants have the same weight, if "SKAT", the rarest variants have the highest weights: weights = -0.4*log10(MAF) )

max.maf.causal

The maf threshold to consider a rare variant (set at 0.01 by default). Only variants with a MAF upper than this threshold will be kept to compute RVAT power. If simus.haplos="liability", variants with a MAF upper this threshold will have a weight of 0

maf.filter

The MAF filter to apply after the simulations to select rare variants to keep for RVAT power analysis. By default corresponds to max.maf.causal

p.causal

The percentage of causal variants, only needed if simus.haplos = "liability"

p.protect

The proportion of protective variants among causal variants, only needed if simus.haplos = "liability"

h2

The variance explained by the gene, only needed if simus.haplos = "liability"

prev

A vector with the prevalence in each group of individuals, only needed if simus.haplos = "liability"

normal.approx

TRUE/FALSE: whether to use the normal approximation to compute thresholds. Set at TRUE by default, only needed if simus.haplos = "liability"

size

The sizes of each group of individuals

verbose

Whether to display information about the function actions

alpha

The significance level to compute the power

RVAT

On which RVAT among 'CAST', 'WSS' and 'SKAT' to compute power (only needed if power.type="simulations"

SKAT.method

Which method to use to compute SKAT ppower, i.e. permutations or theoretical moments (cf SKAT documentation)

simus.haplos

Which method to simulate the data, if simus.haplos="freqs", rbm.haplos.freqs() is used, otherwise rbm.haplos.thresholds() is used.

replicates

The number of simulations to perform to estimate the power

rep.by.causal

The number of time causal variants will be sampled

cores

How many cores to use for moments computation, set at 10 by default

Details

Simulations are perfromed accordingly to rbm.haplos.thresholds() or rbm.haplos.freqs(). Please refer to the corresponding manuals for more details on the simulation procedures. Variants are filtered after the simulations to keep only the rare ones, defined by maf.filter. By defaut, it corresponds to max.maf.causal is used. To disable this filter, set maf.filter at 0.5.

Value

Power values of RVAT

Examples

  #Simulations of 5 groups of individuals with haplotypes frequencies
  #from the 5 EUR populations

  #Load LCT dataset for haplotype matrix
  data(LCT.haplotypes)
  #Haplotypes for the variants in the LCT gene in the EUR population
  LCT.gene.hap <- LCT.hap[which(LCT.sample$super.population=="EUR"),
                          which(LCT.snps$pos>=136545410 & LCT.snps$pos<=136594750)]

  #Individuals from EUR
  LCT.sample.EUR <- subset(LCT.sample, super.population=="EUR")
  #Matrix of haplotypic frequencies
  LCT.freqs <- sapply(unique(LCT.sample.EUR$population), function(z)
                      ifelse(LCT.sample.EUR$population==z,
                             1/table(LCT.sample.EUR$population)[z], 0))

  #Simulation of genetic data for five groups of 50 individuals
  rbm.haplos.power(haplos=LCT.gene.hap, freqs=LCT.freqs, size=rep(50,5), 
                   replicates=5, rep.by.causal = 5,  RVAT = "CAST", 
                   alpha = c(0.001,2.5e-6), cores = 1)