Hi everyone. I'm hoping to compare two sets of p-values distributed across all SNPs in a genome to look for loci where one p-value (or effect size) is more significant than the other. More specifically, we're interested in sites where a p-value from one contrast is not significant, suggesting an interesting lack of differentiation, while another is. I've been able to visualize the results I'm looking for, but I'm struggling to find an appropriate test that's statistically robust.
For background, I am working with a complex population genomics dataset where we tested for rapid adaptation to a toxin over time with 3 treatment groups:
(1) "E" = negative control population not exposed to toxin
(2) "PA" = positive control population pre-adapted to toxin
(3) "S" = naive experimental population exposed to the toxin throughout the experiment
We ran pooled population sequencing to test for changes in allele frequencies over time, and ran many GLMs to test for differences in treatments, time-point, and interactions between the two. For simplicity here, lets just consider data from the final time-point, where we can assume S population have already adapted to the toxin.
Now, we're left with several sets of p-values representing all pairwise treatment contrasts for every site in the genome. The good news is that we see strong signals of adaptation in both the PA and S populations, relative to E negative control populations. The bad news, the adaptation was so strong, it resulted in major selective sweeps that bury many of our target peaks of interest when looking at manhattan plots.
Fortunately, our experiment was designed in such a way that we should be able to compare results across contrasts to ask very specific questions across our dataset. For example, we've run non-parametric Spearman rank correlations and built RRHOs with log-transformed p-values to test for similarities across contrasts. We found that p-values resulting from "E vs S" contrasts are highly correlated with p-values from "E vs PA" contrasts, suggesting some convergent rapid evolution.
Now, we are looking to do the inverse to pull out p-value differences between contrasts to help control for the strong LD resulting from genomic draft.
We want to find sites with differences in p-values or effect sizes between: "PA vs E" and "PA vs S". In other words, we want to find regions where the normalized manhattan peak height in "PA vs E" is significantly higher than peaks from "PA vs S".
This is because S and E populations had the same starting allele frequencies at the beginning of the experiment, and we expect that S allele frequencies moved towards PA populations due to rapid adaptation to the experimental toxin at our sites of interest. PA populations came from the same source population too, but were adapted prior to the start of the experiment, therefore are likely divergent through a combination of selection and drift.
If we are able to run some statistical test, we are hoping it would effectively mask shared variation from E to look for regions where PA and S have become more similar, i.e. have a significant LACK of a peak in "PA vs S" relative to any pairwise contrast including E. By looking for regions where PA and S are most similar relative to E, we should be able to detect sites experiencing convergent evolution.
To this point, I've been able to visualize what we're look for by simply normalizing -log10(p) peak height between manhattans to the maximum peak height in either plot, then subtracting "PA vs S" -log10(p) values from "PA vs E" -log10(p) values.
In the attached plot, you can see the extreme LD in each individual GLM contrast (top 2 panels) that made peak-detection a nightmare, but we're most interested in the peaks that remain after masking by "PA vs S" (bottom panel). You can ignore the red points, but if you're curious those are just some missense SNPs. Can anyone recommend a statistically robust way to test for these differences in p-values/effect sizes beyond just asking which sites remain significant after masking?
Thanks in advance for any help you can provide.
0 answers
No answers yet.
Log in to answer this question.