Parametric mean differences

Author

Gibran Hemani

Published

September 20, 2024

Background

Testing if the mean selection coefficient amongst SNPs is larger than expected, given a background null distribution. How do we get a beta, standard error and p-value for the mean difference?

Generate data

library(dplyr)

Attaching package: 'dplyr'
The following objects are masked from 'package:stats':

    filter, lag
The following objects are masked from 'package:base':

    intersect, setdiff, setequal, union
null_draws <- rnorm(10000, mean=3, sd=8)
true_value <- 12
hist(null_draws, breaks=100)
abline(v=true_value, col="red", lwd=5)

Empirical p-value

empirical_p <- sum(true_value < null_draws) / length(null_draws)
empirical_p
[1] 0.1333

Parametric p-value

null_mean <- mean(null_draws)
null_sd <- sd(null_draws)
parametric_p <- 1 - pnorm(true_value, mean=null_mean, sd=null_sd)
parametric_p
[1] 0.1319797

beta and standard error

b <- true_value - null_mean
z <- qnorm(parametric_p, lower.tail=FALSE)
se <- b/z

tibble(b=b, se=se, parpval=parametric_p, emp_pval=empirical_p)
# A tibble: 1 × 4
      b    se parpval emp_pval
  <dbl> <dbl>   <dbl>    <dbl>
1  8.94  8.00   0.132    0.133
sessionInfo()
R version 4.4.1 (2024-06-14)
Platform: aarch64-apple-darwin20
Running under: macOS Sonoma 14.6.1

Matrix products: default
BLAS:   /Library/Frameworks/R.framework/Versions/4.4-arm64/Resources/lib/libRblas.0.dylib 
LAPACK: /Library/Frameworks/R.framework/Versions/4.4-arm64/Resources/lib/libRlapack.dylib;  LAPACK version 3.12.0

locale:
[1] en_US.UTF-8/en_US.UTF-8/en_US.UTF-8/C/en_US.UTF-8/en_US.UTF-8

time zone: Europe/London
tzcode source: internal

attached base packages:
[1] stats     graphics  grDevices utils     datasets  methods   base     

other attached packages:
[1] dplyr_1.1.4

loaded via a namespace (and not attached):
 [1] digest_0.6.35     utf8_1.2.4        R6_2.5.1          fastmap_1.2.0    
 [5] tidyselect_1.2.1  xfun_0.44         magrittr_2.0.3    glue_1.7.0       
 [9] tibble_3.2.1      knitr_1.47        pkgconfig_2.0.3   htmltools_0.5.8.1
[13] rmarkdown_2.27    generics_0.1.3    lifecycle_1.0.4   cli_3.6.2        
[17] fansi_1.0.6       vctrs_0.6.5       compiler_4.4.1    tools_4.4.1      
[21] pillar_1.9.0      evaluate_0.23     yaml_2.3.8        rlang_1.1.3      
[25] jsonlite_1.8.8    htmlwidgets_1.6.4