# Is it possible to go back and do pairwise comparisons following an LR test, if so what is the correct way

snp1, snp2, and snp3 are "versions" of geneA. I want to determine if geneA significantly effects disease. I have sex as a potentially confounding factor. If geneA affects disease I want to do pairwise comparisons to see which snp affects disease most significantly. I have tens of thousands of these genes to look at, so I am trying to limit type 1 error by only looking at pairwise comparisons in genes that have significant overall effect.

This is my dummy data for geneA

data <- data.frame("snp1"=c(runif(n=150, min=0,max=2),
c(runif(n=50, min=0,max=2))),
"snp2"=c(runif(n=50, min=0, max=.2),
runif(n=50, min=0, max=.2),
runif(n=50, min=1.5, max=2),
runif(n=50, min=1.5, max=2)),
"snp3"=c(runif(n=50, min=0, max=.2),
runif(n=50, min=0, max=.2),
runif(n=50, min=1.5, max=2),
runif(n=50, min=1.5, max=2)),
"sex"=runif(n=50, min=0, max=1),
"disease"=c(rbinom(150, 1, 0.1),
rbinom(50, 1, 0.9)))


This is my likelihood ratio test approach for the whole gene, with snp1 left out as reference:

multi_snp_full <- glm(disease ~ snp2 + snp3 + sex, data=data, family="binomial")
null <- glm(disease ~ sex, data=data, family="binomial")
anova( null, multi_snp_full , test="Chisq")


If I wanted to go back and see how snp2 affect disease would it just be this? (eg no LR test):

single_snp_test <- glm(disease ~ snp2 + sex, data=data, family="binomial")
summary(single_snp_test )

• Since you have snp3 in your full model as well, you should look at the effect of snp2 in multi_snp_full Commented Nov 20, 2020 at 16:29
• >you should look at the effect of snp2 in multi_snp_full I get that snp3 is in the full model, but I don't really get what you mean by the second half of you comment. Could you please clarify what you mean? Commented Nov 20, 2020 at 16:37
• Wait are you saying that if I was intersted in determining the effect of snp2 or snp3 by themselves I would just go: summary(multi_snp_full) and look at the effect and p value for snp2 or snp3 in the readout? So no need to make another model like I did with single_snp_test Commented Nov 20, 2020 at 16:47
• Yes, that is correct. Commented Nov 20, 2020 at 17:43