rvtests icon indicating copy to clipboard operation
rvtests copied to clipboard

Dose the result have a difference between Rvtest and R SKAT package?

Open jeong2624 opened this issue 3 years ago • 0 comments

Hi! I'm interested in rare-variant test, such as burden test, SKAT, SKATO and so on. As I used SKAT test in Rvtest and R SKAT package, The results of p-value are different by using Rvtest and R SKAT package.
I think the results are same, because of the same test. I don't know why the results are different. So, what should I do?

example.vcf data are below:

#CHROM POS ID REF ALT QUAL FILTER INFO FORMAT P1 P2 P3 P4 P5 P6 P7 P8 P9 1 1 . A G 100 . . GT 0/1 0/0 0/0 1/1 0/1 0/0 0/1 0/1 0/0 1 2 . A G 100 . . GT 0/0 0/0 0/0 0/0 0/0 0/0 0/1 0/0 0/1 1 3 . A G 100 . . GT 0/0 0/0 0/0 0/0 0/1 0/0 0/0 0/0 0/0

pheno data are below: fid iid fatid matid sex y1 y2 y3 y4 P1 P1 0 0 0 1.911 -1.465 -0.817 1 P2 P2 0 0 2 2.146 -2.451 -0.178 2 P3 P3 0 0 2 1.086 -1.194 -0.899 1 P4 P4 0 0 2 0.704 -1.052 -0.237 1 P5 P5 0 0 1 2.512 -3.085 -2.579 1 P6 P6 0 0 2 1.283 -1.294 -0.342 1 P7 P7 0 0 1 2.384 1.070 -1.522 2 P8 P8 0 0 1 3.004 0.680 -1.875 1 P9 P9 0 0 1 0.714 -0.116 -1.604 1

setFile data are below: 1:1-3

The my command in Rvtest rvtest --noweb --inVcf ./example.vcf.gz --pheno ./pheno --pheno-name y2 --setFile ./setFile --impute hwe --out ./skat/output --kernel skat

The my coomand in R script library(SKAT) pheno = read.table("/home/pjw/NGS_Project/Rare_excercise/rvtest_result/excercise/pheno",header = T) pheno = pheno[,c("y2")] pheno = as.matrix(pheno) obj = SKAT_Null_Model(pheno ~ 1, out_type="C", n.Resampling = 10000)

library(vcfR) geno = read.vcfR("./example.vcf") gt <- extract.gt(geno) gt <- extract.gt(geno, element = 'GT', as.numeric = TRUE) write.table(gt, "./geno_matrix")

geno = read.table("/home/pjw/NGS_Project/Rare_excercise/rvtest_result/excercise/geno_matrix") geno = as.matrix(t(geno))

SKAT(geno,obj,impute.method = "fixed", weights.beta = c(0.05,1,25))$p.value

jeong2624 avatar Aug 17 '20 07:08 jeong2624