Changes

From Genome Analysis Wiki
Jump to navigationJump to search

SeqShop: Genetic Association Analysis Practical, June 2014 (view source)

Revision as of 12:43, 19 June 2014

3,991 bytes added ,  12:43, 19 June 2014
no edit summary
Line 216: Line 216:  
You can look also visualize the results by QQ-plot and Manhattan plot
 
You can look also visualize the results by QQ-plot and Manhattan plot
    +
<div class="mw-collapsible mw-collapsed" style="width:400px">
 +
View QQ plots
 +
<div class="mw-collapsible-content" style="width:800px">
 +
[[File:Single.epacts.qq.png]]
 +
</div>
 +
</div>
    +
<div class="mw-collapsible mw-collapsed" style="width:400px">
 +
View Manhattan plots
 +
<div class="mw-collapsible-content" style="width:800px">
 +
[[File:Single.epacts.mh.png]]
 +
</div>
 +
</div>
   −
=== Visualizing Zoom Plot ===
+
Also, you can create a zoom plot focusing on the region of interest
   −
* Use ENCODE GM12878 cell type to diplsy variants in the locusZoom-like plot
+
$EPACTS/bin/epacts-zoom --vcf $OUT/assoc/snps.anno.vcf.gz --pos 22:36995620 --prefix $OUT/assoc/single
   −
  $S5/epacts/bin/epacts-zoom --vcf ~/out/snps/chr7.filtered.PASS.beagled.anno.vcf.gz --pos 7:117149147 --prefix ~/out/assoc/single.b.score --cellType Gm12878
+
<div class="mw-collapsible mw-collapsed" style="width:400px">
 +
View Zoom Plots
 +
<div class="mw-collapsible-content" style="width:800px">
 +
[[File:Single.zoom.22.36995620.png]]
 +
</div>
 +
</div>
 +
 
 +
If you want to run EMMAX, you first need to create a kinship matrix
 +
 
 +
  $EPACTS/bin/epacts-make-kin --vcf $OUT/assoc/snps.anno.vcf.gz --min-maf 0.01 --out $OUT/assoc/snps.anno.kinf --run 2 --chr 22
 +
 
 +
And run EMMAX test specifying the kinship matrix
 +
 
 +
$EPACTS/bin/epacts-single --ped $IN/assoc/seqshop.ped --vcf $OUT/assoc/snps.anno.vcf.gz --out $OUT/assoc/emmax --region 22:36000000-37000000 --test q.emmax --pheno PHENO --run 2 --kinf $OUT/assoc/snps.anno.kinf
 +
 
 +
Then the results may look similar to previous ones.
 +
 
 +
cat $OUT/assoc/emmax.epacts.top5000
   −
Copy the output file single.b.score.zoom.7.117149147.pdf using WinSCP or other file transfer software
      
== Run Groupwise Test ==
 
== Run Groupwise Test ==
   −
* Create a set of markers to test
+
To run group-wise test such as gene-level burden test, you need to create a marker group file. If your VCF is already annotated, you can create a group file by running
 +
 
 +
$EPACTS/bin/epacts make-group --vcf $OUT/assoc/snps.anno.vcf.gz --out $OUT/assoc/snps.anno.grp --nonsyn
 +
 
 +
The group file is simply a list of marker per group name, as shown below.
 +
 
 +
cat $OUT/assoc/snps.anno.grp
 +
APOL1 22:36655735_G/A 22:36657740_G/A 22:36661330_G/A 22:36661566_G/A 22:36661646_G/A 22:36661891_G/A 22:36661906_A/G
 +
APOL2 22:36623731_T/C 22:36623920_G/A 22:36629466_T/A 22:36633107_C/A
 +
APOL3 22:36537763_C/T 22:36537798_G/A 22:36556768_G/A 22:36556823_G/T
 +
APOL4 22:36587154_G/T 22:36587202_G/A 22:36587223_G/T 22:36587346_C/T 22:36587511_C/T 22:36587704_T/C 22:36587886_C/T 22:36593714_G/A 22:36597744_A/C 22:36598049_C/G 22:36598058_T/C 22:36598081_A/T
 +
APOL5 22:36122356_G/A 22:36122380_T/A 22:36122930_C/T 22:36123083_C/T 22:36124860_C/G
 +
FOXRED2 22:36900271_T/C 22:36900806_A/G
 +
MYH9 22:36681163_G/C 22:36684354_T/C 22:36710183_T/C
 +
Metazoa_SRP 22:36711990_C/G
 +
RBFOX2 22:36424450_A/C
 +
RP4-633O19__A.1 22:36792162_G/A
 +
 
 +
If you have your own annotation, you can create your own burden test unit by modifying this file.
 +
 
 +
If you want to run a collapsing burden test (CMC), run the following command
 +
 
 +
$EPACTS/bin/epacts group --ped $IN/assoc/seqshop.ped --vcf $OUT/assoc/snps.anno.vcf.gz --out $OUT/assoc/group.collapse --test b.collapse --groupf $OUT/assoc/snps.anno.grp --pheno PHENO --run 2
   −
  $S5/epacts/bin/epacts make-group --vcf ~/out/snps/chr7.filtered.PASS.beagled.anno.vcf.gz --out ~/out/snps/chr7.filtered.PASS.beagled.anno.grp --nonsyn
+
You can view the results by examining the output file
 +
  cat $OUT/assoc/group.collapse.epacts
   −
* Run SKAT-O test
+
<div class="mw-collapsible mw-collapsed" style="width:400px">
 +
View Output file
 +
<div class="mw-collapsible-content" style="width:800px">
 +
#CHROM BEGIN END MARKER_ID NS FRAC_WITH_RARE NUM_ALL_VARS NUM_PASS_VARS NUM_SING_VARS PVALUE STATRHO
 +
22 36655735 36661906 22:36655735-36661906_APOL1 62 0.14516 7 4 0 0.42748 1
 +
22 36623731 36633107 22:36623731-36633107_APOL2 62 0.064516 4 1 0 0.038657 NA
 +
22 36537763 36556823 22:36537763-36556823_APOL3 62 0.080645 4 2 0 0.40634 0
 +
22 36587154 36598081 22:36587154-36598081_APOL4 62 0.14516 12 4 0 0.67891 0
 +
22 36122356 36124860 22:36122356-36124860_APOL5 62 0.1129 5 2 0 0.15055 0.3
 +
22 36900271 36900806 22:36900271-36900806_FOXRED2 NA NA 2 0 0 NA NA
 +
22 36681163 36710183 22:36681163-36710183_MYH9 62 0.032258 3 1 0 1 NA
 +
22 36711990 36711990 22:36711990-36711990_Metazoa_SRP NA NA 1 0 0 NA NA
 +
22 36424450 36424450 22:36424450-36424450_RBFOX2 62 0.032258 1 1 0 1 NA
 +
22 36792162 36792162 22:36792162-36792162_RP4-633O19__A.1 NA NA 1 0 0 NA NA
 +
</div>
 +
</div>
   −
$S5/epacts/bin/epacts group --ped $S5/examples/index/chr7.CFTR.ped --vcf ~/out/snps/chr7.filtered.PASS.beagled.anno.vcf.gz --out ~/out/assoc/group.skat.o --groupf ~/out/snps/chr7.filtered.PASS.beagled.anno.grp --test skat --skat-o --run 2
+
You can run SKAT-O test in a similar way, but with a special tag
   −
  cat ~/out/assoc/group.skat.o.epacts | cut -f 4,6,10,11
+
  $EPACTS/bin/epacts group --ped $IN/assoc/seqshop.ped --vcf $OUT/assoc/snps.anno.vcf.gz --out $OUT/assoc/group.skato --test skat --skat-o --groupf $OUT/assoc/snps.anno.grp --pheno PHENO --run 2
   −
* Run Variable Threshold test
+
And view output files
−
  $S5/epacts/bin/epacts group --ped $S5/examples/index/chr7.CFTR.ped --vcf ~/out/snps/chr7.filtered.PASS.beagled.anno.vcf.gz --out ~/out/assoc/group.VT --groupf ~/out/snps/chr7.filtered.PASS.beagled.anno.grp --test VT --run 2
+
   
 +
cat $OUT/assoc/group.skato.epacts
   −
  cat ~/out/assoc/group.VT.epacts | cut -f 4,6,10,13
+
<div class="mw-collapsible mw-collapsed" style="width:400px">
 +
View Output file
 +
<div class="mw-collapsible-content" style="width:800px">
 +
#CHROM BEGIN END MARKER_ID NS FRAC_WITH_RARE NUM_ALL_VARS NUM_PASS_VARS NUM_SING_VARS PVALUE STATRHO
 +
22 36655735 36661906 22:36655735-36661906_APOL1 62 0.14516 7 4 0 0.42748 1
 +
22 36623731 36633107 22:36623731-36633107_APOL2 62 0.064516 4 1 0 0.038657 NA
 +
  22 36537763 36556823 22:36537763-36556823_APOL3 62 0.080645 4 2 0 0.40634 0
 +
22 36587154 36598081 22:36587154-36598081_APOL4 62 0.14516 12 4 0 0.67891 0
 +
22 36122356 36124860 22:36122356-36124860_APOL5 62 0.1129 5 2 0 0.15055 0.3
 +
22 36900271 36900806 22:36900271-36900806_FOXRED2 NA NA 2 0 0 NA NA
 +
22 36681163 36710183 22:36681163-36710183_MYH9 62 0.032258 3 1 0 1 NA
 +
22 36711990 36711990 22:36711990-36711990_Metazoa_SRP NA NA 1 0 0 NA NA
 +
22 36424450 36424450 22:36424450-36424450_RBFOX2 62 0.032258 1 1 0 1 NA
 +
22 36792162 36792162 22:36792162-36792162_RP4-633O19__A.1 NA NA 1 0 0 NA NA
 +
</div>
 +
</div>
Retrieved from "http://genome.sph.umich.edu/wiki/Special:MobileDiff/11006"

Navigation menu