1
$\begingroup$

I have next SNPs data.frame structure

> str(SNPs)
'data.frame':   1703 obs. of  4 variables:
 $ group: Factor w/ 2 levels "A","B": 1 1 1 1 1 1 1 1 1 1 ...
 $ rs1  : Factor w/ 3 levels "D/D","I/D","I/I": 1 1 2 3 3 2 1 1 1 1 ...
 $ rs2  : Factor w/ 3 levels "a/a","a/b","b/b": 3 3 2 3 3 2 2 3 3 2 ...
 $ rs3  : Factor w/ 3 levels "G/G","G/T","T/T": 2 1 2 1 1 3 1 1 2 1 ...
 ...other rs


> head(SNPs)
  group rs1 rs2 rs3 ...other rs
1     A D/D b/b G/T
2     A D/D b/b G/G
3     A I/D a/b G/T
4     A I/I b/b G/G
5     A I/I b/b G/G
6     A I/D a/b T/T

For example, I noticed that rs5 4b/4a and rs6 G/G in group A very often occur together (see below). In group B they very often occur together too. So I want to know - is it statistical regularity or not.

I can create table with all pairs in both groups

> SNPs$rs5_rs6 <- paste(SNPs$rs5, SNPs$rs6)
> tmp <- table(SNPs$rs5_rs6, SNPs$group)
> tmp

              A   B
  4a/4a G/G   1  20
  4b/4a G/G  31  83
  4b/4a G/T  14  51
  4b/4a T/T   1   0
  4b/4b G/G  37 106
  4b/4b G/T  35 119
  4b/4b T/T  11  31

So, now I need compare (find p-value) group A and group B: 4a/4a G/G in group A vs 4a/4a G/G in group B, 4b/4a G/G in group A vs 4b/4a G/G in group B, 4b/4a G/T in group A vs 4b/4a G/T in group B and etc. Only this pairs/

How can I do that? I need p-values for all

> tmp

              A   B
  4a/4a G/G   1  20 - p-value?
  4b/4a G/G  31  83 - p-value?
  4b/4a G/T  14  51 - p-value?
  4b/4a T/T   1   0 - p-value?
  4b/4b G/G  37 106 - p-value?
  4b/4b G/T  35 119 - p-value?
  4b/4b T/T  11  31 - p-value?
$\endgroup$
5
  • $\begingroup$ If you can restate your question in terms of explicit null and alternative hypotheses it will be much easier for people to help you determine an appropriate test to find p-values $\endgroup$ Jul 23, 2018 at 14:41
  • $\begingroup$ @heathobrien chisq - it is my choice. Problem is not to choice test. Problem is to calculate it for pairs in tmp $\endgroup$
    – n.osennij
    Jul 23, 2018 at 15:04
  • $\begingroup$ If your question is how to apply a chisq test to the above table, the answer is chisq.test(tmp$A, tmp$B) $\endgroup$ Jul 23, 2018 at 15:22
  • $\begingroup$ heathobrien, my question is how to apply chisq test to tmp and get 7 p-values: compare4a/4a G/G in group A vs 4a/4a G/G in group B and get first p-value, compare 4b/4a G/G in group A vs 4b/4a G/G in group B and get the second p-value and etc. $\endgroup$
    – n.osennij
    Jul 23, 2018 at 15:26
  • 1
    $\begingroup$ OK. Now we're getting somewhere. See my answer $\endgroup$ Jul 23, 2018 at 16:01

1 Answer 1

2
$\begingroup$

You can get standardised residuals from chisq.test and convert these to p-values, but it's not recommended. The frequencies of each genotype depend on the counts of all other genotypes, so the p-values aren't independent and they aren't corrected for multiple testing, but if you must:

pnorm(-abs(chisq.test(tmp[,c('A','B')])$stdres)))
$\endgroup$
8
  • $\begingroup$ hm... What other options are possible? whether the differences between the two groups were statistically significant for all pairs (rs5+rs6 and all pairs - 4a/4a G/G, 4b/4a G/G and etc)? $\endgroup$
    – n.osennij
    Jul 23, 2018 at 16:07
  • $\begingroup$ If I were you, I would just focus on the overall chisq p-value (corrected for the number of SNP pairs tested). The frequency of allele 2 depends on the frequency of allele 1 and the frequency of the genotypes depends on the allele frequencies, so I don't think it makes much sense to test each genotype independently $\endgroup$ Jul 23, 2018 at 16:18
  • 1
    $\begingroup$ Note that all you're doing is looking at differences in linkage disequilibrium between groups. There are a lot of tools out there to that efficiently $\endgroup$ Jul 23, 2018 at 16:41
  • 1
    $\begingroup$ Fishers exact test doesn't give you residuals. Having values less than 5 is probably the least of your worries if you're trying to obtain p-values for each combination of genotypes. $\endgroup$ Jul 23, 2018 at 20:36
  • 1
    $\begingroup$ Let us continue this discussion in chat. $\endgroup$ Jul 23, 2018 at 20:45

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge that you have read and understand our privacy policy and code of conduct.

Not the answer you're looking for? Browse other questions tagged or ask your own question.