Athurva Gore: Difference between revisions

From ZhangLabWiki
Jump to navigation Jump to search
>Ajgore
>Ajgore
Line 324: Line 324:
  codeStart will be set to the starting point of coding, codeEnd will be set to the end point of coding.
  codeStart will be set to the starting point of coding, codeEnd will be set to the end point of coding.
**This should allow coding and non-coding regions to be identified by the probe generation script.
**This should allow coding and non-coding regions to be identified by the probe generation script.
**Finished script; ran script on all RefSeq exons!
**All refseq exons are now present in a mysql database on laptop.
*'''Test with probe generation code'''
**Using "SELECT INTO FILE," a set of exons was created to test the probe generation code.
**Uploaded file to pgp.oxf, testing it now.
***Due to "brute-force" algorithm for genes, it seems to be taking a long time; make sure script is working.

Revision as of 21:12, 12 September 2008

Lab Rotation

  • Will add daily entries once schedule settles down a bit.

August

August 4th - August 10th

  • Received CpG island probe code from Dr. Zhang.
  • Need to make modifications discussed in meeting:
  • Necessary modifications:
    • Remove bisulfite conversion function - DONE
      • Create dummy function that simply does not do a replacement, but creates fwd and rev strands
      • Commented code to find/replace
      • Question: parts we are looking at; do they have N's in them? Probably shouldn't.
    • Remove CG score penalty - DONE
      • in function getProbeList
      • remove section starting with "my $CG_penalty"
      • set CG_penalty as always zero for quick fix, remove CG finding code later
    • Add function to check sequence uniqueness in human genome (build hash table ahead of time, check Dr. Zhang's functions)
      • Looks like function oligoFreq does this; currently just uses log of each side
      • Modify probe score calculation to favor segment B being more unique
      • See function getProbeList
    • Modify score calculation so that Tm of segment A should be high and Tm of segment B should be low
      • At the moment, Tm score is used in p_score calculation
      • Simply modify to adjust so that A has high Tm, B has low Tm
      • See function getProbeList
      • Specific thresholds?
    • Exon database
      • Ask about our current UCSC data; it comes in format on PGP server, or is this modified?
      • Found Exon-Intron Database; uses GenBank
      • Show Dr. Zhang the readme

August 11th - August 17th

  • Will copy from notebook

August 18th - August 24th

  • Will copy from Notebook

August 27th

  • Have code to read and analyze Cosmic set data; file is "readCosmic.pl"
  • Comparison of Cosmic set data to Dr. Li's data (Harvard Collaborator)
    • Dr. Li's data showed that the last base pair of the ligation arm was important.
    • Average expression was calculated using the "makeHist.pl" file.
      • This file averages the abundance of the probes that used each base.
    • The results were slightly different.
      • For the COSMIC set, it was found that for the proximal ligation arm, G was preferred, followed by C, A, and T.
      • This makes intuitive sense, as the bond between G and C is stronger than the bond between A and T.
      • The results were also much less pronounced than Dr. Li's results; the difference in base performance is much smaller.
Cosmic Set Results
Base: A T G C
Proximal Extension Arm: 98.2276826484018 105.4814504701790 97.8751164168995 105.1120477176380
Distal Extension Arm: 104.9745183977170 107.4126226040130 106.5526921141540 88.8525735059634
Proximal Ligation Arm: 95.1374247734938 83.7531848897751 135.360909000299 107.876484056988
Distal Ligation Arm: 103.1547820103830 113.0664396095970 96.5690024678530 96.7110820814283
  • File:ProxDist Base Cosmic.png
    • Analyze these results later; perhaps add in a preference on the Proximal Ligation Arm to favor G or C over A or T in the score.
  • Implementation of new probe set algorithm.
    • Algorithm is summarized above.
    • Implement in Perl using hash table in getCpGIslands.pl file.
  • Use R to analyze Cosmic Data Set
    • Downloaded and installed GNU R and R-commander (GUI)
      • Can load CSV and tab-separated files.
    • lm command can do a linear model fit to the data.
      • May need to take the log of several pieces of data however in order to improve fit (if data is exponentially correlated).
      • Initial lm results for Cosmic set show a very poor correlation...probably need to take log of several pieces of data set.
    • Earlier, Dr. Zhang's analysis on CES22K set used the following factors:
      • ProbeDist: distribution of oligo abundance determined by Solexa sequencing. - Don't have this for the Cosmic set
      • H1_Tm, H2_Tm: H1 (ligation arm), H2 (extension arm), 12bp sub-sequences of H1 and H2. - Normalized for the Cosmic set; very close to 62 C in all cases
      • H1R_dist, H2L_dist: K-mer distribution of H1&H2 12bp sub-sequences in the human genome. - Have this
      • Core_25bp_genomic_homology: weighted score of the 25bp homologous sequences in the human genome, mapped by SOAP with <=3bp mismatches. - Not sure what this is...ask Dr. Zhang

August 28th

  • TO DO TODAY:
    • Take log of various factors in R, see if this affects lm results.
    • See if any other combinations of factors produce better fits in R.
    • Implement algorithm to select probes in Perl.
  • R Analysis of Cosmic Data Set
    • Tried taking the log of all factors aside from Tm, checked fit.
    • Fit is still not very good
      • Very low R^2 value
      • Checked by plotting in Octave; fit does not match up to data well at all.
    • Perhaps a different, non-linear, non-logarithmic model is needed?
    • One issue with the Cosmic Data Set is that there is no data for oligonucleotide bias on the chip.
      • Agilent's chips have some inherent bias in what molecules are present and not present.
      • This bias was measured for the CES22k set that Dr. Zhang worked on but not for the Cosmic set.
      • This may be causing issues with fitting a curve to the Cosmic set, as this bias cannot be used in the fit.
    • Came up with some ideas while working from home.
      • Perhaps instead of log or linear fit, use a normal distribution type fit?
      • Could fit a normal distribution to each factor and try and standardize it; this might result in a better overall fit.
      • I am currently looking into how to easily do this in R or Octave; shouldn't be too hard to do a fit this way.
    • Received data for CES22K set
      • A preliminary look at the data in Octave makes it seem that R^2 values will remain very low.
      • Confirmed that fit is not that great in Octave.
      • However, data for bias is present; this means that this can be factored into the analysis.
      • Work on this tomorrow.

August 29th

  • TO DO TODAY:
    • Make some slides for Cosmic and CES22k sets to illustrate findings for linear and log fits
    • Hopefully find a good way to fit normal distributions to the data; create some plots for this
    • Test CES22k set for proximal and distal bases as well (just as with Dr. Li's data and Cosmic set above)
  • Analysis of Cosmic and CES22k Data Set
    • R^2 value of fit is very low due to large number of data points.
      • Idea from Dr. Zhang; check high and low end points for input factors (maybe outside 1 std. deviation?)
      • This might increase Signal-to-Noise ratio for the data.
    • Was doing some reading on R linear model modules; stumbled across effects module.
      • effects uses results from an lm command and generates "effect plots"
      • It generates the effects of a set of variables on an output variable in individual graph form
      • Red lines in graphs represent a 95% confidence interval
      • A paper on the effects module can be found here: http://www.jstatsoft.org/v08/i15/paper
      • Need to investigate this module and make sure that any assumptions it makes are valid for the data we are using.
        • It seemed that it was being used in a variety of applications online that used linear models.
    • Used the effects module on Dr. Zhang's earlier CES22k fit to determine effects of various factors on expression output.
      • The results for this are very interesting, and provide a much nicer picture than just the fit data.



  • Proximal and Distal Arm Base Comparison
    • The above analysis (using makeHist_CES.pl) was repeated for the CES22k set.
    • The results are below.
CES22k Set Results
Base: A T G C
Proximal Extension Arm: 305.94 326.29 236.49 247.1
Distal Extension Arm: 283.08 290.35 263.52 257.27
Proximal Ligation Arm: 244.76 163.81 402.33 253.13
Distal Ligation Arm: 280.15 303.46 245.07 263.64
  • File:ProxDist Base CES22k.png
    • Similar to Cosmic set results; G and C are preferred at the Proximal Ligation site over A and T.
    • Differences are more pronounced, like Dr. Li's data.
    • Some differences are present in the other bases as well...possibly analyze for significance (in Dr. Li's data, Cosmic data, and this data)?
  • Slides to illustrate findings so far
    • Will make tonight (8-29) and email out.

September

September 1

  • Made notes for powerpoint slides at home; however, came up with new things to do before powerpoint is finished.
  • Will be presenting on Thursday.

September 2

  • TO DO:
    • Add in GC content, AT content for arms and target region for effects plots.
    • Make slides; show each factor's effect graphs with explanation and a comparison of Proximal/Distal Bases (including Dr. Li's data)
      • Describe how these curves will be used to affect scoring in next set of probes
      • Take how score is currently calculated into account, suggest modifications?
      • Decide how best to scale everything.
  • PowerPoint
    • Analyzed data and made preliminary powerpoint of results.
    • Found that:
      • 1) For the proximal ligation character, G is preferred, followed by C, A, and then T.
      • 2) A high Tm for the ligation arm is good.
      • 3) A low Tm for the extension arm is good.
      • 4) Smaller targets are captured with higher efficiency.
      • 5) Low GC content is better for both arms according to the Cosmic set. I still need to look at the CES22K set.
      • 6) Conclusions cannot really be drawn from the frequency data for the Cosmic set or the CES22K set, as the confidence interval is very very large. In general, the only thing I can see is that the ligation arm frequency should be higher.
    • Here is the finished PowerPoint: Media:Cosmic_CES22K_Presentation.ppt
    • Still need to add target GC content data and arm GC content from CES22K set.
  • Incorporating Effect Plots into Scoring
    • Effect plots have allowed each factor to have an equation fitting its value to output.
    • This allows these equations to be directly used in determining probe design.
    • Perhaps compare current scoring method (that uses match between oligos A and B in code) to output as well?
    • Then, would be able to determine how much to scale this factor as well as the other, already fit factors.

September 3

  • TO DO:
    • Analyze arm GC content from CES22K Set and add to Powerpoint.
    • Add arm length data from CES22K Set and Cosmic Set to Powerpoint.
    • Implement new scoring system.
  • Data Analysis and PowerPoint Presentation
    • Analyzed arm length from Cosmic set; found that longer arms (on both sides) result in a better capture.
      • Will add this to PowerPoint and to scoring.
    • Need to extract and analyze GC Content from CES22K Set.
      • Code to do this is running right now. After this, graphs will be created in R and added to PowerPoint and scoring system.
    • Turns out arm lengths were normalized to 20 for this set; cannot derive anything from this because of this.
    • GC results are different; it was found that the ligation arm having a high GC raised expression, while the extension arm having a high GC hurt expression.
      • In the Cosmic set, higher GC in either arm led to decreased expression.
      • I am not sure why this is the case; look into how probes were constructed for tips.
  • New Scoring System
    • For scoring, intercept of fit does not really matter, since all scores are just going to be compared to each other. We therefore just need to apply the slope.
    • New Score:
      • p_score = m_LTm * Left_Tm + m_RTm * Right_Tm + m_LUnique * L_Unique + m_RUnique * R_Unique + m_TLength * T_Length + m_LLength * L_Length + m_RLength * R_Length + m_LGC * GC_Content_Left + m_RGC * GC_Content_Right + proxLigFactor
    • This scoring should produce much better probes than previous iterations, as all of these factors will be fitted.
    • Additional factors can then be analyzed as well.
    • Added new scoring system to GetCpGProbes.pl file. Just need to decide exactly which slopes to use for probe scoring.
      • Will probably average results from Cosmic and CES22K (after making sure that other factors didn't cause differences; if this is the case, a weighted average will be used).

September 4

  • TO DO:
    • Give presentation on data analysis to Ida, Alice, and Dr. Zhang.
    • Finish up code for new probe design code (add in slopes from R analysis, finalize algorithm to select all probe sets).
  • Presentation
    • Added some finishing touches to PowerPoint; will go over it to make sure that I know it well.
    • Make sure that it looks good; make any changes necessary for style.
  • New Probe Design Code
    • Added in framework for scoring yesterday.
    • Just need to add in slopes today.
    • Maybe have the code read them from a file? Won't have to edit code every time, can instead just edit a "scoring method" file.
    • Will average slopes from Cosmic and CES22K where both are present.
    • Will use average offset as well for the ligation base character.

September 5

  • TO DO:
    • Plug in final choices for slopes into Perl code.
    • Make sure new probe selection algorithm works.
    • Make sure new scoring works.
    • Try code out on some sequences.
  • Slope Choices
    • Based on the data in the above PowerPoint presentation, the following slopes were chosen for the code.
$m_LLength=18.51978; #Ligation arm length - COSMIC
$m_RLength=13.74286; #Extension arm length - COSMIC
$m_TLength=-1.89955; #Target length - COSMIC
$m_LTm=8.3619; #Ligation Tm - CES22K
$m_RTm=-5.3790; #Extension Tm - CES22K
$m_LGC=-326.17; #Ligation GC content - COSMIC - Maybe reduce a little if this is too much?
$m_RGC=-217.60; #Extension GC content - COSMIC - Maybe reduce a little if this is too much?
$m_LUnique=0; #Ignored for now.
$m_RUnique=0; #Ignored for now.
$g_Offset=268.845; #Offset for Proximal Ligation Base 'G'- Average of COSMIC and CES22K
$a_Offset=169.95; #Offset for Proximal Ligation Base 'A'- Average of COSMIC and CES22K
$t_Offset=123.78; #Offset for Proximal Ligation Base 'T'- Average of COSMIC and CES22K
$c_Offset=180.505; #Offset for Proximal Ligation Base 'C'- Average of COSMIC and CES22K
    • This data was added to the GetCpgProbes.pl file.
  • Test probe selection algorithm and new scoring system
    • Chromosome 20, being the smallest chromosome, was used for testing (speed purposes).
    • A section was chosen for selection.
      • Upload test file, link here
    • The code was tested and debugged.

September 6-10

  • TO DO:
    • Finish Perl probe generation script
  • Worked on Perl script
  • Learned several issues with how Perl passes arrays.
    • Implemented workarounds in code; used split and join in order to pass multiple arrays as multiple scalars.
  • Algorithm to make sure all possible probe combinations were used was implemented using recursion
    • Several bugs in recursion were present, causing infinite loops.
    • Fixed this; recursive algorithm to generate probe sets now works.
  • Looked around for list of Cosmic set genes
    • Found in Cosmic_3760_exons.list.txt
      • Located in /home/kun/MIP/exonome on pgp.oxf freelogy server
    • Format of set is:
exonID		chromosome	start		end
NM_032211-E2	chr10		100000811	100000923
NM_032211-E3	chr10		100001312	100001449
NM_032211-E4	chr10		100002099	100002215
...
...

September 11

  • TO DO:
    • Modify Perl probe generation script to read files in Cosmic set format (like above)
    • Post new Perl script here
    • Run tests, compare new results to old Cosmic set results
    • Show Dr. Zhang the code
  • Cosmic File Reading
    • Need to modified Perl script so that Cosmic set type files can be read in successfully.
    • Also need to make sure that Gene and Exon IDs are output by the final program; otherwise it will be very difficult to keep track of probes.
    • These modifications have been made; an in-file setting now controls what file format to load from.
  • Fixing errors with some probes
    • It was discovered that some probe sets were causing errors while others were not; this was related to the way Perl's "return" function works. The probeList variables were not undefined, as they contained a single number with the error code for not finding any compatible probes.
    • This was fixed.
  • New Perl Script
  • NEW: Redundant Sequences
    • Now, need to figure out a better way to deal with redundant portions of various genes.
    • Basically, different genes that include the same exon.
    • New data structure?
      • Perhaps create exonIDs for various genes? Then, list genes as a combination of exons in this database?
      • Prevents overlap.
    • Added refGene data from UCSC to a mysql table called refGene on laptop.
      • Data can be ordered by chromosome, txStart and end regions.
    • However, multiple exons are stored in single rows.
    • A script could be used (PHP?) to extract each exon and create a new database table.
    • Redundant Sequence Script Idea:
      • Use PHP.
      • Scan through table; grab all values from exonStarts and exonEnds (these are our exons.)
      • For each exon, also grab the "name" field, the "chr" field, and the "strand" field.
      • Reprint out new exons using a basic incrementing ID. Also reprint out gene data in a new table as well.
Table Exon - Fields:
exonID (autoincrement)     chr     strand     start     end
Table ExonLink - Fields:
geneID     exonID     position
NOTE: There will be multiple entries for each geneID.
      In order to figure out which exons to sequence for a given gene, all exonIDs will be pulled out of the database.
      Probes for each exonID will be designed.
      Example:  SELECT exonID FROM TABLE gene WHERE geneID=geneid_value
  • After this is finished, probes can be ordered using new code.

September 12

  • TO DO:
    • Finish script that creates new tables.
  • Script to create refExon and refExonLink
    • Created a PHP script that can build the two new tables (refExon and refExonLink).
    • In order to differentiate between coding and non-coding regions, in the refExonLink table, I have added two columns: codeStart and codeEnd.
    • Basically, these two values allow us to figure out what is coding and what is not in the following way.
    • Cases:
      • Out of bounds
xxxxxx |   |      OR      |   | xxxxxx
Both codeStart and codeEnd will be set to -1.
      • In bounds
| xxxxxx |
Both codeStart and codeEnd will be set to -10.
      • Overlap
xxxx|xxx   |
codeStart will be set to the starting point of coding, codeEnd will be set to -10.

|   xxxx|xxx
codeStart will be set to -10, codeEnd will be set to the end point of coding.

xx|xxx|xx
codeStart will be set to the starting point of coding, codeEnd will be set to the end point of coding.
    • This should allow coding and non-coding regions to be identified by the probe generation script.
    • Finished script; ran script on all RefSeq exons!
    • All refseq exons are now present in a mysql database on laptop.
  • Test with probe generation code
    • Using "SELECT INTO FILE," a set of exons was created to test the probe generation code.
    • Uploaded file to pgp.oxf, testing it now.
      • Due to "brute-force" algorithm for genes, it seems to be taking a long time; make sure script is working.