Changes

From Genome Analysis Wiki
Jump to navigationJump to search

SeqShop: Ancestry On Your Own Genome, May 2015 (view source)

Revision as of 14:57, 21 May 2015

641 bytes removed ,  14:57, 21 May 2015
→‎Visualizing Ancestry
Line 15: Line 15:       −
<div class="mw-collapsible mw-collapsed" style="width:700px">
  −
'''Only if your snpcalling is still running, let Mary Kate know, expand, and follow these instructions'''
  −
<div class="mw-collapsible-content">
  −
* Detach from screen
  −
Ctrl-a d
  −
* Start new screen
  −
screen
  −
* Reset environment variables (replace SampleXX with your sample name)
  −
export SAMPLE=SampleXX
  −
source /net/seqshop-server/home/mktrost/seqshop/setupSS.txt
  −
* When you detach from screen and reattach, you will have multiple screen sessions
  −
** The ancestry one will have today's date, the snpcall one will have an earlier date.
  −
*** Mary Kate can explain how to reattach to a specific session when we get there.
  −
</div>
  −
</div>
  −
  −
  −
'''Instructions for everyone: Setup LASER'''
     −
Source the LASER setup:
+
'''Setup LASER'''
−
  source /net/seqshop-server/home/chaolong/LASER-Tutorial/setup.txt
+
  export SSREF=/net/seqshop-server/home/mktrost/seqshop/singleSample/ref
−
(SAMPLE & OUT were previously set in your screen session.)
+
export LASER=~/seqshop/output/ancestry/LASER-2.01
 +
(SAMPLE were previously set in your screen session.)
    
After setting this, also do
 
After setting this, also do
Line 48: Line 31:  
=== Step 1: bam --> pileup ===
 
=== Step 1: bam --> pileup ===
   −
  $GC/bin/samtools mpileup -q 30 -Q 20 -f $REF/human.g1k.v37.fa -l $HGDP/HGDP_938.bed $OUT/bams/${SAMPLE}.recal.bam > $OUT/ancestry/${SAMPLE}.HGDP.pileup
+
  $GC/bin/samtools mpileup -q 30 -Q 20 -f $SSREF/gotcloud.ref/hs37d5.fa -l $SSREF/HGDP/HGDP_938.bed $OUT/bams/${SAMPLE}.recal.bam > $OUT/ancestry/${SAMPLE}.HGDP.pileup
    
This step takes ~1 hour for a genome sequenced at 17X.
 
This step takes ~1 hour for a genome sequenced at 17X.
Line 55: Line 38:  
  Ctrl-a d
 
  Ctrl-a d
   −
While that's running, let's go run the Association Analysis Practical: [[SeqShop: Genetic Association Analysis Practical, December 2014]]
+
While that's running, let's go look at Calling results: [[ SeqShop: Calling Your Own Genome, May 2015]]
    
=== Checking if Pileup finished ===
 
=== Checking if Pileup finished ===
Line 63: Line 46:       −
If not, let's go to Indel: [[SeqShop:_Calling_Your_Own_Genome,_December_2014#Reviewing_Indel_Results]]
+
If not, let's go to Indel: [[SeqShop:_Calling_Your_Own_Genome,_May_2015#Reviewing_Indel_Results]]
       
If yes, let's continue below:
 
If yes, let's continue below:
 +
 
=== Step 2: pileup --> seq ===
 
=== Step 2: pileup --> seq ===
    
  python $LASER/pileup2seq/pileup2seq.py \
 
  python $LASER/pileup2seq/pileup2seq.py \
−
  -m $HGDP/HGDP_938.site \
+
  -m $SSREF/HGDP/HGDP_938.site \
 
  -o $OUT/ancestry/$SAMPLE.HGDP \
 
  -o $OUT/ancestry/$SAMPLE.HGDP \
 
  $OUT/ancestry/$SAMPLE.HGDP.pileup
 
  $OUT/ancestry/$SAMPLE.HGDP.pileup
Line 79: Line 63:  
This step will take about 5-6 minutes.
 
This step will take about 5-6 minutes.
   −
  $LASER/laser -g $HGDP/HGDP_938.geno -c $HGDP/HGDP_938.RefPC.coord -s $OUT/ancestry/$SAMPLE.HGDP.seq -K 20 -k 4 -M 0.8 -o $OUT/ancestry/$SAMPLE.HGDP &
+
  $LASER/laser -g $SSREF/HGDP/HGDP_938.geno -c $SSREF/HGDP/HGDP_938.RefPC.coord -s $OUT/ancestry/$SAMPLE.HGDP.seq -K 20 -k 4 -M 0.8 -o $OUT/ancestry/$SAMPLE.HGDP &
    
View the results:
 
View the results:
Line 92: Line 76:     
Generate the plot:
 
Generate the plot:
−
  Rscript plotHGDP.r $HGDP/HGDP_938.RefPC.coord $OUT/ancestry/${SAMPLE}.HGDP.SeqPC.coord  
+
  Rscript plotHGDP.r $SSREF/HGDP/HGDP_938.RefPC.coord $OUT/ancestry/${SAMPLE}.HGDP.SeqPC.coord  
    
Take a look:
 
Take a look:
Line 118: Line 102:     
For European ancestry:
 
For European ancestry:
−
  $LASER/laser -g $HGDP/HGDP.633K.euro.geno -c $HGDP/HGDP.633K.euro.RefPC.coord -s $OUT/ancestry/$SAMPLE.HGDP.seq -K 20 -k 4 -M 0.8 -o $OUT/ancestry/$SAMPLE.HGDP-Euro &
+
  $LASER/laser -g $SSREF/HGDP/HGDP.633K.euro.geno -c $SSREF/HGDP/HGDP.633K.euro.RefPC.coord -s $OUT/ancestry/$SAMPLE.HGDP.seq -K 20 -k 4 -M 0.8 -o $OUT/ancestry/$SAMPLE.HGDP-Euro &
    
For East Asian ancestry:
 
For East Asian ancestry:
−
  $LASER/laser -g $HGDP/HGDP.633K.easia.2.geno -c $HGDP/HGDP.633K.easia.2.RefPC.coord -s $OUT/ancestry/$SAMPLE.HGDP.seq -K 20 -k 4 -M 0.8 -o $OUT/ancestry/$SAMPLE.HGDP-EAsia &
+
  $LASER/laser -g $SSREF/HGDP/HGDP.633K.easia.2.geno -c $SSREF/HGDP/HGDP.633K.easia.2.RefPC.coord -s $OUT/ancestry/$SAMPLE.HGDP.seq -K 20 -k 4 -M 0.8 -o $OUT/ancestry/$SAMPLE.HGDP-EAsia &
    
For Central/South Asian ancestry:
 
For Central/South Asian ancestry:
−
  $LASER/laser -g $HGDP/HGDP.633K.csasia.2.geno -c $HGDP/HGDP.633K.csasia.2.RefPC.coord -s $OUT/ancestry/$SAMPLE.HGDP.seq -K 20 -k 4 -M 0.8 -o $OUT/ancestry/$SAMPLE.HGDP-CSAsia &
+
  $LASER/laser -g $SSREF/HGDP/HGDP.633K.csasia.2.geno -c $SSREF/HGDP/HGDP.633K.csasia.2.RefPC.coord -s $OUT/ancestry/$SAMPLE.HGDP.seq -K 20 -k 4 -M 0.8 -o $OUT/ancestry/$SAMPLE.HGDP-CSAsia &
      Line 144: Line 128:     
Generate the plot for European ancestry:
 
Generate the plot for European ancestry:
−
  Rscript plotHGDP.r $HGDP/HGDP.633K.euro.RefPC.coord $OUT/ancestry/${SAMPLE}.HGDP-Euro.SeqPC.coord  
+
  Rscript plotHGDP.r $SSREF/HGDP/HGDP.633K.euro.RefPC.coord $OUT/ancestry/${SAMPLE}.HGDP-Euro.SeqPC.coord  
 
or for East Asian ancestry
 
or for East Asian ancestry
−
  Rscript plotHGDP.r $HGDP/HGDP.633K.easia.2.RefPC.coord $OUT/ancestry/${SAMPLE}.HGDP-EAsia.SeqPC.coord  
+
  Rscript plotHGDP.r $SSREF/HGDP/HGDP.633K.easia.2.RefPC.coord $OUT/ancestry/${SAMPLE}.HGDP-EAsia.SeqPC.coord  
 
or for Central/South Asian ancestry
 
or for Central/South Asian ancestry
−
  Rscript plotHGDP.r $HGDP/HGDP.633K.csasia.2.RefPC.coord $OUT/ancestry/${SAMPLE}.HGDP-CSAsia.SeqPC.coord  
+
  Rscript plotHGDP.r $SSREF/HGDP/HGDP.633K.csasia.2.RefPC.coord $OUT/ancestry/${SAMPLE}.HGDP-CSAsia.SeqPC.coord  
    
Take a look:
 
Take a look:
Line 154: Line 138:     
== Return to Calling your own Genome ==
 
== Return to Calling your own Genome ==
−
Go to [[SeqShop:_Calling_Your_Own_Genome,_December_2014#Reviewing_Indel_Results]] and pick up where we left off.
+
Go to [[SeqShop:_Calling_Your_Own_Genome,_May_2015#Reviewing_Indel_Results]] and pick up where we left off.
Retrieved from "http://genome.sph.umich.edu/wiki/Special:MobileDiff/13284...13410"

Navigation menu