rvtest

rvtest: Rare Variants tests

下载

git clone https://github.com/zhanxw/rvtests

单变量测试

rvtest --inVcf input.vcf --pheno phenotype.ped --out output --single wald,score

# input. vcf
#CHROM  POS     ID      REF     ALT     QUAL    FILTER  INFO    FORMAT  P1      P2      P3      P4      P5 
1       1       .       A       G       100     .       .       GT      0/1     0/0     0/0     1/1     0/1
1       2       .       A       G       100     .       .       GT      0/0     0/0     0/0     0/0     0/0 

# phenotype.ped
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

/Users/chengkai/Documents/05_software/rvtests/rvtests_macos/executable/rvtest --inVcf 01new_cohort_ID.vcf --pheno pheno1.txt --out output --single wald,score

image.png

分组测试

• Burden tests: group variants, which are usually less than 1% or 5% rare variants, for association tests. The category includes: CMC test, Zeggini test, Madsen-Browning test, CMAT test, and rare-cover test.
• Variable threshold tests: group variants under different frequency thresholds.
• Kernel methods: suitable to tests rare variants having different directions of effects. These includes SKAT test and KBAC test.

rvtest --inVcf input.vcf --pheno phenotype.ped --out output --geneFile refFlat_hg19.txt.gz --burden cmc --vt price --kernel skat,kbac

# vcf 需要修改，去掉前面的chr
def change_cohort_ID():
    with open("/Users/chengkai/Documents/lishuyuan/code/new_method/filter/cohort/rvtest/cohort_ID.vcf") as f:
        for i in f.readlines():
            if i.startswith("#"):
                print(i,end="")
            else:
                i=i.strip().split()
                i[0]=i[0][3:]
                print("\t".join(i))
change_cohort_ID()

# vcf 需要建立索引文件 
bgzip vcf
tabix -p vcf

/Users/chengkai/Documents/05_software/rvtests/rvtests_macos/executable/rvtest --inVcf 01new_cohort_ID.vcf --pheno pheno1.txt --out output --geneFile refFlat_hg19.txt.gz --burden cmc --vt price --kernel skat,kbac,skato

image.png

• 输出参数详解

Gene    RANGE   N_INFORMATIVE   NumVar  NumPolyVar      Q       rho     Pvalue
OR4F5   1:69090-70008   125     2       2       107.481 1       0.693376
SAMD11  1:861120-879961 125     36      35      579.602 1       0.8867
NOC2L   1:879582-894679 125     48      47      3354.42 1       1

• NumVar: number of variants in the gene(or site)

• NumPolyVar: number of Polymorphic Genotypes

• NonRefSite: non-reference site, genotype is not exactly 0.0

• 问题与讨论

image.png