Matt:LabNotes/2013-8-20
Jump to navigation
Jump to search
CA12k Capture Analysis (MiSeq_130325)[edit]
- Here is improved analysis of the MiSeq_130325 data
- Previously used Bowtie2 custom parameters, which now I know are not optimal
- Previously aligned to hg19 and now aligned to probelist used to order the probes
- The old analysis can be found in the commented out section at the bottom of the "edit" page or Matt:LabNotes/MiSeq_130325_Analysis
Building bowtie2 index for CA12k probes[edit]
- First converted probelist to fasta file with H1, H2,and 100bp linking region rearranged using CAProbelist_to_fasta.pl: Feb2013_200bp_probes_to_order.txt -> CAprobes_to_order_H1H2.fa
- Fasta reference made to match reads: RevComp(H1) + RevComp(H2) + RevComp(Barcode)
- Build index: genome-miner:$ bowtie2-build CAprobes_to_order_H1H2.fa CAprobes_H1H2
Processing reads[edit]
- Shorten Read1 reads to 146bp to match reference length with h1h2_first146_only.pl (nothing fancy, just substr($_,0,146) for base call and base quality strings)
D1-cDNA-RNaseA_S2_L001_R1_001.fastq -> cDNA-RNaseA_R1_146shortened.fastq D1-cDNAwRNaseA_S3_L001_R1_001.fastq -> cDNAwRNaseA_R1_146shortened.fastq D1-gDNA_S1_L001_R1_001.fastq -> gDNA_R1_146shortened.fastq
Map reads with default Bowtie2[edit]
bowtie2 --phred33 -x CAprobes_H1H2 -q gDNA_R1_146shortened.fastq > gDNA_R1_146shortened_CAprobesH1H2.sam & 4476041 reads; of these: 4476041 (100.00%) were unpaired; of these: 129983 (2.90%) aligned 0 times 4322633 (96.57%) aligned exactly 1 time 23425 (0.52%) aligned >1 times 97.10% overall alignment rate
bowtie2 --phred33 -x CAprobes_H1H2 -q cDNA-RNaseA_R1_146shortened.fastq > cDNA-RNaseA_R1_146shortened_CAprobesH1H2.sam 2> cDNA-RNase_stderr.txt & 2852904 reads; of these: 2852904 (100.00%) were unpaired; of these: 81850 (2.87%) aligned 0 times 2740392 (96.06%) aligned exactly 1 time 30662 (1.07%) aligned >1 times 97.13% overall alignment rate
bowtie2 --phred33 -x CAprobes_H1H2 -q cDNAwRNaseA_R1_146shortened.fastq > cDNAwRNaseA_R1_146shortened_CAprobesH1H2.sam 2> cDNAwRNase_stderr.txt & 2091417 reads; of these: 2091417 (100.00%) were unpaired; of these: 55021 (2.63%) aligned 0 times 2010743 (96.14%) aligned exactly 1 time 25653 (1.23%) aligned >1 times 97.37% overall alignment rate
Sort and filter (MAPQ>=8) aligned reads[edit]
- MAPQ=8 corresponds to 16% chance alignment is wrong, but manual checking showed that alignments with MAPQ=8 matched well
samtools view -bS gDNA_R1_146shortened_CAprobesH1H2.sam | samtools sort - gDNA_R1_146shortened_CAprobesH1H2_sorted samtools view -h -F 4 -q 8 gDNA_R1_146shortened_CAprobesH1H2_sorted.bam > gDNA_R1_146shortened_CAprobesH1H2_sorted_filtered.sam
samtools view -bS cDNA-RNaseA_R1_146shortened_CAprobesH1H2.sam | samtools sort - cDNA-RNaseA_R1_146shortened_CAprobesH1H2_sorted samtools view -h -F 4 -q 8 cDNA-RNaseA_R1_146shortened_CAprobesH1H2_sorted.bam > cDNA-RNaseA_R1_146shortened_CAprobesH1H2_sorted_filtered.sam
samtools view -bS cDNAwRNaseA_R1_146shortened_CAprobesH1H2.sam | samtools sort - cDNAwRNaseA_R1_146shortened_CAprobesH1H2_sorted samtools view -h -F 4 -q 8 cDNAwRNaseA_R1_146shortened_CAprobesH1H2_sorted.bam > cDNAwRNaseA_R1_146shortened_CAprobesH1H2_sorted_filtered.sam
Count reads aligned for each probe[edit]
- Used Count_reads_from_bowtie2_SAM.pl to count number of appearances of each probe in column 3 (Reference Name) of sam file
- CountofcDNA-RNaseA_CAprobesH1H2_filtered.txt -> total count for all probes: 2,554,856
- CountofcDNAwRNaseA_CAprobesH1H2_filtered.txt -> total count for all probes: 1,881,008
- CountofgDNA_CAprobesH1H2_filtered.txt -> total count for all probes: 3,976,822
Histograms and Results[edit]
- Only used cDNA-RNaseA since previous regression analysis showed with or without RNaseA were essentially the same
- "Good" probes were the set of probes that had a capture efficiency within one order of magnitude
- "Brain optimized" probes were the set of probes that captured
File:CA12k gDNA Capture.JPGFile:CA12k cDNA-RNaseA Capture.JPG
- 3,141 Unique "good" genes
- 2,615 Unique "brain optimized" genes
- 4,185 "Brain optimized" and "good" probes
- 1,963 "Brain optimized" and "good" genes
- 1,849 of which are brain genes
- 189 of which are stem cell genes
Details of genes/probes can be found in CA12k_Capture_Analysis.xlsx
Comparison of CA12k to Agi26k Results Summary[edit]
' | CA12k | Agi26k0gap | Agi26k20gap |
Good probes | 7,900 | 6,596 | 5,403 |
Good genes | 3,141 | 2,628 | 2,529 |
Brain-optimized probes | 6,412 | 8,640 | 8,438 |
Brain-optimized genes | 2,615 | 3,140 | 3,133 |
Good & Brain-optimized probes | 4,185 | 4,822 | 4,006 |
Good & Brain-optimized genes | 1,963 | 2,075 | 2,075 |
Good & Brain-optimized brain genes | 1,849 | 2,002 | 1,996 |
Good & Brain-optimized stem cell genes | 189 | 138 | 151 |