s2p: s2p

Description Usage Arguments Value Examples

View source: R/genotype_sim.R

Description

Inverse function for p2s. Converts expected values of SD(Z) (ie, tau, sigma_1, or sigma_2) into underlying 'population' odds-ratio distribution. Assume that, across some set of SNPs in, population log-odds ratios between two phenotypes A and B are normally distributed with standard deviation p. If a GWAS is performed between a group of samples with phenotype A of size n1 and a group of samples of phenotype B of size n2, and Z-scores calculated for each SNP with standard deviation s, this function recovers p.

Usage

1
s2p(s, n1, n2, mafs = NULL, NS = 50000)

Arguments

s

standard deviation of underlying log odds ratio distribution

n1

number of samples in group 1

n2

number of samples in group 2

mafs

use to specify distribution of average MAF (across both groups). Default is MAF~U(0.01,0.5)

NS

number of SNPs to use in simulation

Value

expected value of observed standard deviation of Z scores.

Examples

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
14
15
16
17
18
19
p_underlying=0.15; # underlying standard deviation of odds-ratio distribution
nsnp=10000; # number of SNPs
n1=500; n2=1000 # number of samples in group 1, group 2

ors=exp(rnorm(nsnp,mean=0,sd=p_underlying)); # (simulated) underlying odds ratios

pmaf_1=runif(nsnp,0.05,0.5); # (simulated) underlying population minor allele frequencies in group 1
pmaf_2=ors*pmaf_1/(1-pmaf_1+(ors*pmaf_1)) # underlying population minor allele frequencies in group 2
plot(pmaf_2*(1-pmaf_1)/(pmaf_1*(1-pmaf_2)),ors); abline(0,1,col="red") # check

omaf_1=(rbinom(nsnp,n1,p=pmaf_1)+rbinom(nsnp,n1,p=pmaf_1))/(2*n1) # (simulated) observed diploid minor allele frequencies in group 1
omaf_2=(rbinom(nsnp,n2,p=pmaf_2)+rbinom(nsnp,n2,p=pmaf_2))/(2*n2) # (simulated) observed diploid minor allele frequencies in group 1

test=apply(cbind(omaf_1,omaf_2),1,function(x) prop.test(2*c(n1,n2)*x,2*c(n1,n2))$p.value) # chi-square tests between frequencies under H0: pmaf_1==pmaf_2

obs_z=-qnorm(test/2)*sign(omaf_2-omaf_1) # observed Z scores

p_underlying # observed SD of Z scores
s2p(sd(obs_z),n1,n2) # predicted

jamesliley/subtest documentation built on May 18, 2019, 11:21 a.m.