<?xml version='1.0'?><rss version="2.0" xmlns:dc="http://purl.org/dc/elements/1.1/" xmlns:georss="http://www.georss.org/georss" xmlns:atom="http://www.w3.org/2005/Atom" >
<channel>
	<title><![CDATA[BOL: All site blogs]]></title>
	<link>https://bioinformaticsonline.com/blog/all?offset=150</link>
	<atom:link href="https://bioinformaticsonline.com/blog/all?offset=150" rel="self" type="application/rss+xml" />
	<description><![CDATA[]]></description>
	
	<item>
	<guid isPermaLink="true">https://bioinformaticsonline.com/blog/view/36425/check-product-if-already-exist-in-to-cart-using-woocommerce</guid>
	<pubDate>Sun, 29 Apr 2018 00:27:22 -0500</pubDate>
	<link>https://bioinformaticsonline.com/blog/view/36425/check-product-if-already-exist-in-to-cart-using-woocommerce</link>
	<title><![CDATA[Check product if already exist in to cart using woocommerce.]]></title>
	<description><![CDATA[<p>Check product if already exist in to cart using woocommerce.<br /><br />add_action( 'woocommerce_add_to_cart_validation', 'check_product_added_to_cart', 10, 3 );<br /><br />function check_product_added_to_cart( $passed, $product_id, $quantity) {<br /><br />foreach (WC()-&gt;cart-&gt;get_cart() as $cart_key =&gt; $cart_item ){<br />// if products are already in cart:<br />if( $cart_item['product_id'] == $product_id ) {<br />// Set the verification variable to "not passed" (false)<br />$passed = false;<br />// (Optionally) Displays a notice if product(s) are already in cart<br />wc_add_notice( '&lt;strong&gt;' . $btn['label'] . '&lt;/strong&gt; ' . __( 'This product is already in your cart.', 'woocommerce' ), 'error' );<br />// Stop the function returning "false", so the products will not be added again<br />return $passed;<br />}<br />}<br />return $passed;<br />}</p><p>More at https://phpajaxhtml.blogspot.in/2017/11/check-product-if-already-exist-in-to.html</p>]]></description>
	<dc:creator>Ram Yash Pal</dc:creator>
</item>
<item>
	<guid isPermaLink="true">https://bioinformaticsonline.com/blog/view/36413/insert-data-through-ajax-into-mysql-database</guid>
	<pubDate>Wed, 25 Apr 2018 12:43:51 -0500</pubDate>
	<link>https://bioinformaticsonline.com/blog/view/36413/insert-data-through-ajax-into-mysql-database</link>
	<title><![CDATA[Insert data through ajax into MySql database.]]></title>
	<description><![CDATA[<p>Insert data through ajax into MySql database.1. Create form.php and copy below code into file.</p><blockquote><p>&lt;script type="text/javascript"&gt;<br /> <br /> $(document).ready(function(){<br /> $('#submit').click(function(){<br /> var name= $('#name').val();<br /> var email= $('#email').val();<br /> var sdatatring='name1='+ name +'&amp;email1='+ email;<br /> $.ajax({<br /> type:"POST",<br /> url:"insert.php",<br /> data:sdatatring,<br /> cache: false,<br /> success:function(result){<br /> alert(result);</p><p>}});</p><p>});</p><p>});<br />&lt;/script&gt;<br /> &lt;form method="post" action="" name="frm"&gt;<br /> Name:&lt;input type="text" name="name" id="name" value=""&gt;&lt;br&gt;<br /> Email:&lt;input type="text" name="email" id="email" value=""&gt;&lt;br&gt;<br /> &lt;input type="button" name="submt" id="submit" value="submit" /&gt;<br /> <br /> &lt;/form&gt;</p><p>2. Create insert.php and copy below code into file.<br /> &lt;?php<br />print_r($_POST);<br />$con=mysql_connect("localhost","root","");<br />mysql_select_db('dbname');<br />mysql_query('insert into tablename('colname')values(value)');<br /> ?&gt;<br />Note:You need to include jQuery library in form.php.</p></blockquote><p>More at</p><p>https://phpajaxhtml.blogspot.in/2018/03/insert-data-through-ajax-into-mysql.html</p>]]></description>
	<dc:creator>Ram Yash Pal</dc:creator>
</item>
<item>
	<guid isPermaLink="true">https://bioinformaticsonline.com/blog/view/36411/create-custom-widget-area-in-wordpress</guid>
	<pubDate>Wed, 25 Apr 2018 12:35:53 -0500</pubDate>
	<link>https://bioinformaticsonline.com/blog/view/36411/create-custom-widget-area-in-wordpress</link>
	<title><![CDATA[Create custom widget area in wordpress.]]></title>
	<description><![CDATA[<p>https://phpajaxhtml.blogspot.in/2018/04/create-custom-widget-in-wordpress.html</p>]]></description>
	<dc:creator>Ram Yash Pal</dc:creator>
</item>
<item>
	<guid isPermaLink="true">https://bioinformaticsonline.com/blog/view/36211/project-based-approach-to-improve-bioinformatics-education-with-skilled-and-meaningful-access-to-omics-data</guid>
	<pubDate>Wed, 11 Apr 2018 13:31:42 -0500</pubDate>
	<link>https://bioinformaticsonline.com/blog/view/36211/project-based-approach-to-improve-bioinformatics-education-with-skilled-and-meaningful-access-to-omics-data</link>
	<title><![CDATA[Project-based approach to improve bioinformatics education with skilled and meaningful access to omics data]]></title>
	<description><![CDATA[<p>Pine Biotech has been collaborating with Loyola University of New Orleans on piloting a new approach to bioinformatics education using the intuitive and logic-drive bioinformatics platform T-BioInfo.</p><p>https://edu.t-bio.info/collaborative-model-bioinformatics-education-combining-biologically-inspired-bioinformatics-project-based-learning/</p>]]></description>
	<dc:creator>eliabrodsky</dc:creator>
</item>
<item>
	<guid isPermaLink="true">https://bioinformaticsonline.com/blog/view/35923/basic-command-line-to-run-blast</guid>
	<pubDate>Wed, 14 Mar 2018 05:10:34 -0500</pubDate>
	<link>https://bioinformaticsonline.com/blog/view/35923/basic-command-line-to-run-blast</link>
	<title><![CDATA[Basic command-line to run BLAST]]></title>
	<description><![CDATA[<p>&nbsp;</p><p>The goal of this tutorial is to run you through a demonstration of the command line, which you may not have seen or used much before.</p><p>All of the commands below can copy/pasted.</p><div id="install-software"><h2>Install software<a href="http://angus.readthedocs.io/en/2016/running-command-line-blast.html#install-software" title="Permalink to this headline"></a></h2><p>Copy and paste the following commands</p><div><div><pre>sudo apt-get update &amp;&amp; sudo apt-get -y install python ncbi-blast+
</pre></div></div><p>This updates the software list and installs the Python programming language and NCBI BLAST+.</p></div><div id="get-data"><h2>Get Data<a href="http://angus.readthedocs.io/en/2016/running-command-line-blast.html#get-data" title="Permalink to this headline"></a></h2><p>Grab some data to play with. Grab some cow and human RefSeq proteins:</p><div><div><pre>wget ftp://ftp.ncbi.nih.gov/refseq/B_taurus/mRNA_Prot/cow.1.protein.faa.gz
wget ftp://ftp.ncbi.nih.gov/refseq/H_sapiens/mRNA_Prot/human.1.protein.faa.gz
</pre></div></div><p>This is only the first part of the human and cow protein files - there are 24 files total for human.</p><p>The database files are both gzipped, so lets unzip them</p><div><div><pre>gunzip *gz
ls
</pre></div></div><p>Take a look at the head of each file:</p><div><div><pre>head cow.1.protein.faa
head human.1.protein.faa
</pre></div></div><p>These are protein sequences in FASTA format. FASTA format is something many of you have probably seen in one form or another &ndash; it&rsquo;s pretty ubiquitous. It&rsquo;s just a text file, containing records; each record starts with a line beginning with a &lsquo;&gt;&rsquo;, and then contains one or more lines of sequence text.</p><p>Note that the files are in fasta format, even though they end if &rdquo;.faa&rdquo; instead of the usual &rdquo;.fasta&rdquo;. This NCBI&rsquo;s way of denoting that this is a fasta file with amino acids instead of nucleotides.</p><p>How many sequences are in each one?</p><div><div><pre>grep -c '^&gt;' cow.1.protein.faa
grep -c '^&gt;' human.1.protein.faa
</pre></div></div><p>This grep command uses the c flag, which reports a count of lines with match to the pattern. In this case, the pattern is a regular expression, meaning match only lines that begin with a &gt;.</p><p>This is a bit too big, lets take a smaller set for practice. Lets take the first two sequences of the cow proteins, which we can see are on the first 6 lines</p><div><div><pre>head -6 cow.1.protein.faa &gt; cow.small.faa
</pre></div></div></div><div id="blast"><h2>BLAST<a href="http://angus.readthedocs.io/en/2016/running-command-line-blast.html#blast" title="Permalink to this headline"></a></h2><p>Now we can blast these two cow sequences against the set of human sequences. First, we need to tell blast about our database. BLAST needs to do some pre-work on the database file prior to searching. This helps to make the software work a lot faster. Because you installed your own version of the sotware, you need to tell the shell where the software is located. Use the full path and the makeblastdb command:</p><div><div><pre>makeblastdb -in human.1.protein.faa -dbtype prot
ls
</pre></div></div><p>Note that this makes a lot of extra files, with the same name as the database plus new extensions (.pin, .psq, etc). To make blast work, these files, called index files, must be in the same directory as the fasta file.</p><p><br /> blastp [-h] [-help] [-import_search_strategy filename]<br /> [-export_search_strategy filename] [-task task_name] [-db database_name]<br /> [-dbsize num_letters] [-gilist filename] [-seqidlist filename]<br /> [-negative_gilist filename] [-negative_seqidlist filename]<br /> [-entrez_query entrez_query] [-db_soft_mask filtering_algorithm]<br /> [-db_hard_mask filtering_algorithm] [-subject subject_input_file]<br /> [-subject_loc range] [-query input_file] [-out output_file]<br /> [-evalue evalue] [-word_size int_value] [-gapopen open_penalty]<br /> [-gapextend extend_penalty] [-qcov_hsp_perc float_value]<br /> [-max_hsps int_value] [-xdrop_ungap float_value] [-xdrop_gap float_value]<br /> [-xdrop_gap_final float_value] [-searchsp int_value]<br /> [-sum_stats bool_value] [-seg SEG_options] [-soft_masking soft_masking]<br /> [-matrix matrix_name] [-threshold float_value] [-culling_limit int_value]<br /> [-best_hit_overhang float_value] [-best_hit_score_edge float_value]<br /> [-window_size int_value] [-lcase_masking] [-query_loc range]<br /> [-parse_deflines] [-outfmt format] [-show_gis]<br /> [-num_descriptions int_value] [-num_alignments int_value]<br /> [-line_length line_length] [-html] [-max_target_seqs num_sequences]<br /> [-num_threads int_value] [-ungapped] [-remote] [-comp_based_stats compo]<br /> [-use_sw_tback] [-version]</p><p>Now we can run the blast job. We will use blastp, which is appropriate for protein to protein comparisons.</p><div><div><pre>blastp -query cow.small.faa -db human.1.protein.faa
</pre></div></div><p>This gives us a lot of information on the terminal screen. But this is difficult to save and use later - Blast also gives the option of saving the text to a file.</p><div><div><pre>    blastp -query cow.small.faa -db human.1.protein.faa -out cow_vs_human_blast_results.txt
ls
</pre></div></div><p>Take a look at the results using less. Note that there can be more than one match between the query and the same subject. These are referred to as high-scoring segment pairs (HSPs).</p><div><div><pre>less cow_vs_human_blast_results.txt
</pre></div></div><p>So how do you know about all the options, such as the flag to create an output file? Lets also take a look at the help pages. Unfortunately there are no man pages (those are usually reserved for shell commands, but some software authors will provide them as well), but there is a text help output</p><div><div><pre>blastp -help
</pre></div></div><p>To scroll through slowly</p><div><div><pre>blastp -help | less
</pre></div></div><p>To quit the less screen, press the q key.</p><p>Parameters of interest include the -evalue (Default is 10?!?) and the -outfmt</p><p>Lets filter for more statistically significant matches with a different output format:</p><div><div><pre>blastp \
-query cow.small.faa \
-db human.1.protein.faa \
-out cow_vs_human_blast_results.tab \
-evalue 1e-5 \
-outfmt 7
</pre></div></div><p>I broke the long single command into many lines with by &ldquo;escaping&rdquo; the newline. That forward slash tells the command line &ldquo;Wait, I&rsquo;m not done yet!&rdquo;. So it waits for the next line of the command before executing.</p><p>Check out the results with less.</p><p>Lets try a medium sized data set next</p><div><div><pre>head -199 cow.1.protein.faa &gt; cow.medium.faa
</pre></div></div><p>What size is this db?</p><div><div><pre>grep -c '^&gt;' cow.medium.faa
</pre></div></div><p>Lets run the blast again, but this time lets return only the best hit for each query.</p><div><div><pre>blastp \
-query cow.medium.faa \
-db human.1.protein.faa \
-out cow_vs_human_blast_results.tab \
-evalue 1e-5 \
-outfmt 6 \
-max_target_seqs 1
</pre></div></div></div><div id="summary"><h2>Summary<a href="http://angus.readthedocs.io/en/2016/running-command-line-blast.html#summary" title="Permalink to this headline"></a></h2><p>Review:</p><ul>
<li>command line programs such as blast use flags to get information about how and what to do</li>
<li>blast options can be found by typing&nbsp;<cite>blastp -help</cite></li>
<li>break a command up over many lines by using&nbsp;<a href="http://angus.readthedocs.io/en/2016/running-command-line-blast.html#id1">`</a>` to &ldquo;escape&rdquo; the new line</li>
</ul><p>&nbsp;</p><p>Blastn</p><p>blastn [-h] [-help] [-import_search_strategy filename]<br /> [-export_search_strategy filename] [-task task_name] [-db database_name]<br /> [-dbsize num_letters] [-gilist filename] [-seqidlist filename]<br /> [-negative_gilist filename] [-negative_seqidlist filename]<br /> [-entrez_query entrez_query] [-db_soft_mask filtering_algorithm]<br /> [-db_hard_mask filtering_algorithm] [-subject subject_input_file]<br /> [-subject_loc range] [-query input_file] [-out output_file]<br /> [-evalue evalue] [-word_size int_value] [-gapopen open_penalty]<br /> [-gapextend extend_penalty] [-perc_identity float_value]<br /> [-qcov_hsp_perc float_value] [-max_hsps int_value]<br /> [-xdrop_ungap float_value] [-xdrop_gap float_value]<br /> [-xdrop_gap_final float_value] [-searchsp int_value]<br /> [-sum_stats bool_value] [-penalty penalty] [-reward reward] [-no_greedy]<br /> [-min_raw_gapped_score int_value] [-template_type type]<br /> [-template_length int_value] [-dust DUST_options]<br /> [-filtering_db filtering_database]<br /> [-window_masker_taxid window_masker_taxid]<br /> [-window_masker_db window_masker_db] [-soft_masking soft_masking]<br /> [-ungapped] [-culling_limit int_value] [-best_hit_overhang float_value]<br /> [-best_hit_score_edge float_value] [-window_size int_value]<br /> [-off_diagonal_range int_value] [-use_index boolean] [-index_name string]<br /> [-lcase_masking] [-query_loc range] [-strand strand] [-parse_deflines]<br /> [-outfmt format] [-show_gis] [-num_descriptions int_value]<br /> [-num_alignments int_value] [-line_length line_length] [-html]<br /> [-max_target_seqs num_sequences] [-num_threads int_value] [-remote]<br /> [-version]</p><p>DESCRIPTION<br /> Nucleotide-Nucleotide BLAST 2.7.0+</p></div>]]></description>
	<dc:creator>Shruti Paniwala</dc:creator>
</item>
<item>
	<guid isPermaLink="true">https://bioinformaticsonline.com/blog/view/35621/bbtools-for-bioinformatician</guid>
	<pubDate>Thu, 15 Feb 2018 16:45:52 -0600</pubDate>
	<link>https://bioinformaticsonline.com/blog/view/35621/bbtools-for-bioinformatician</link>
	<title><![CDATA[BBTools for bioinformatician !]]></title>
	<description><![CDATA[<p><span></span><br /><strong>BBMap.sh</strong><br /><br /></p><ul>
<li><strong>Mapping Nanopore reads</strong></li>
</ul><p><br /><span>BBMap.sh has a length cap of 6kbp. Reads longer than this will be broken into 6kbp pieces and mapped independently.</span><br /><br /></p><div><div>Code:</div><pre dir="ltr">$ mapPacBio.sh -Xmx20g k=7 in=reads.fastq ref=reference.fa maxlen=1000 minlen=200 idtag ow int=f qin=33 out=mapped1.sam minratio=0.15 ignorequality slow ordered maxindel1=40 maxindel2=400</pre></div><p><br /><span>The "maxlen" flag shreds them to a max length of 1000; you can set that up to 6000. But I found 1000 gave a higher mapping rate.&nbsp;&nbsp;</span><br /><br /></p><ul>
<li><strong>Using Paired-end and single-end reads at the same time</strong></li>
</ul><p><br /><span>BBMap itself can only run single-ended or paired-ended in a single run, but it has a wrapper that can accomplish it, like this:</span><br /><br /></p><div><div>Code:</div><pre dir="ltr">$ bbwrap.sh in1=read1.fq,singletons.fq in2=read2.fq,null out=mapped.sam append</pre></div><p><span>This will write all the reads to the same output file but only print the headers once. I have not tried that for bam output, only sam output</span><br /><br /><span>Note about alignment stats: For paired reads, you can find the total percent mapped by adding the read 1 percent (where it says "mapped: N%") and read 2 percent, then dividing by 2. The different columns tell you the count/percent of each event. Considering the cigar strings from alignment, "Match Rate" is the number of symbols indicating a reference match (=) and error rate is the number indicating substitution, insertion, or deletion (X, I, D).</span><br /><br /></p><ul>
<li><strong>Exact matches when mapping small reads (e.g. miRNA)</strong></li>
</ul><p><br /><span>When mapping small RNA's with BBMap use the following flags to report only perfect matches.</span><br /><br /></p><div><div>Code:</div><pre dir="ltr">ambig=all vslow perfectmode maxsites=1000</pre></div><p><span>It should be very fast in that mode (despite the vslow flag). Vslow mainly removes masking of low-complexity repetitive kmers, which is not usually a problem but can be with extremely short sequences like microRNAs.</span></p><ul>
<li><strong>Important note about BBMap alignments</strong></li>
</ul><p><br /><span>BBMap is always nondeterministic when run in paired-end mode with multiple threads, because the insert-size average is calculated on a per-thread basis, which affects mapping; and which reads are assigned to which thread is nondeterministic. The only way to avoid that would be to restrict it to a single thread (threads=1), or map the reads as single-ended and then fix pairing afterward:</span><br /><br /></p><div><div>Code:</div><pre dir="ltr">bbmap.sh in=reads.fq outu=unmapped.fq int=f
repair.sh in=unmapped.fq out=paired.fq fint outs=singletons.fq</pre></div><p><span>In this case you'd want to only keep the paired output.&nbsp;</span><br /><br /><span>BBSplit is based on BBMap, so it is also nondeterministic in paired mode with multiple threads. BBDuk and Seal (which can be used similarly to BBSplit) are always deterministic.&nbsp;</span><br /><br /><span>--------------------------------------------------------</span><br /><br /><strong>Reformat.sh</strong></p><ul>
<li><strong>Count k-mers/find unknown primers</strong></li>
</ul><p>&nbsp;</p><div><div>Code:</div><pre dir="ltr">$ reformat.sh in=reads.fq out=trimmed.fq ftr=19</pre></div><p><span>This will trim all but the first 20 bases (all bases after position 19, zero-based).</span><br /><br /></p><div><div>Code:</div><pre dir="ltr">$ kmercountexact.sh in=trimmed.fq out=counts.txt fastadump=f mincount=10 k=20 rcomp=f</pre></div><p><span>This will generate a file containing the counts of all 20-mers that occurred at least 10 times, in a 2-column format that is easy to sort in Excel.&nbsp;</span><br /><br /></p><div><div>Code:</div><pre dir="ltr">ACCGTTACCGTTACCGTTAC	100
AAATTTTTTTCCCCCCCCCC	85</pre></div><p><span>...etc. If the primers are 20bp long, they should be pretty obvious.&nbsp;&nbsp;</span></p><ul>
<li><strong>Convert SAM format from 1.4 to 1.3 (required for many programs)</strong></li>
</ul><p>&nbsp;</p><div><div>Code:</div><pre dir="ltr">$ reformat.sh in=reads.sam out=out.sam sam=1.3</pre></div><ul>
<li><strong>Removing N basecalls</strong></li>
</ul><p><br /><span>You can use BBDuk or Reformat with "qtrim=rl trimq=1". That will only trim trailing and leading bases with Q-score below 1, which means Q0, which means N (in either fasta or fastq format). The BBMap package automatically changes q-scores of Ns that are above 0 to 0 and called bases with q-scores below 2 to 2, since occasionally some Illumina software versions produces odd things like a handful of Q0 called bases or Ns with Q&gt;0, neither of which make any sense in the Phred scale.</span></p><ul>
<li><strong>Sampling reads</strong></li>
</ul><p>&nbsp;</p><div><div>Code:</div><pre dir="ltr">$ reformat.sh in=reads.fq out=sampled.fq sample=3000</pre></div><div><div>Code:</div><pre dir="ltr">To sample 10% of the reads:
reformat.sh in1=reads1.fq in2=reads2.fq out1=sampled1.fq out2=sampled2.fq samplerate=0.1

or more concisely:
reformat.sh in=reads#.fq out=sampled#.fq samplerate=0.1

and for exact sampling:
reformat.sh in=reads#.fq out=sampled#.fq samplereadstarget=100k</pre></div><ul>
<li><strong>Changing fasta headers</strong></li>
</ul><p><br /><span>Remove anything after the first space in fasta header.&nbsp;</span><br /><br /></p><div><div>Code:</div><pre dir="ltr"> reformat.sh in=sequences.fasta out=renamed.fasta trd</pre></div><p><span>"trd" stands for "trim read description" and will truncate everything after the first whitespace.</span></p><ul>
<li><strong>Extract reads from a sam file</strong></li>
</ul><p>&nbsp;</p><div><div>Code:</div><pre dir="ltr">$ reformat.sh in=reads.sam out=reads.fastq</pre></div><ul>
<li><strong>Verify pairing and optionally de-interleave the reads</strong></li>
</ul><p>&nbsp;</p><div><div>Code:</div><pre dir="ltr">$ reformat.sh in=reads.fastq verifypairing</pre></div><ul>
<li><strong>Verify pairing if the reads are in separate files</strong></li>
</ul><p>&nbsp;</p><div><div>Code:</div><pre dir="ltr">$ reformat.sh in1=r1.fq in2=r2.fq vpair</pre></div><p><span>If that completes successfully and says the reads were correctly paired, then you can simply de-interleave reads into two files like this:</span><br /><br /></p><div><div>Code:</div><pre dir="ltr">$ reformat.sh in=reads.fastq out1=r1.fastq out2=r2.fastq</pre></div><ul>
<li><strong>Base quality histograms</strong></li>
</ul><p>&nbsp;</p><div><div>Code:</div><pre dir="ltr">$ reformat.sh in=reads.fq qchist=qchist.txt</pre></div><p><span>That stands for "quality count histogram".&nbsp;</span></p><ul>
<li><strong>Filter SAM/BAM file by read length</strong></li>
</ul><p>&nbsp;</p><div><div>Code:</div><pre dir="ltr">$ reformat.sh in=x.sam out=y.sam minlength=50 maxlength=200</pre></div><ul>
<li><strong>Filter SAM/BAM file to detect/filter spliced reads</strong></li>
</ul><p>&nbsp;</p><div><div>Code:</div><pre dir="ltr">$ reformat.sh in=mapped.bam out=filtered.bam maxdellen=50</pre></div><p><span>You can set "maxdellen" to whatever length deletion event you consider the minimum to signify splicing, which depends on the organism.</span><br /><span>-------------------------------------------------------------</span><br /><strong>Repair.sh</strong></p><ul>
<li><strong>"Re-pair" out-of-order reads from paired-end data files</strong></li>
</ul><p>&nbsp;</p><div><div>Code:</div><pre dir="ltr">$ repair.sh in1=r1.fq.gz in2=r2.fq.gz out1=fixed1.fq.gz out2=fixed2.fq.gz outsingle=singletons.fq.gz</pre></div><p><span>--------------------------------------------------------------</span><br /><strong>BBMerge.sh</strong><br /><br /><span>BBMerge now has a new flag - "outa" or "outadapter". This allows you to automatically detect the adapter sequence of reads with short insert sizes, in case you don't know what adapters were used. It works like this:</span><br /><br /></p><div><div>Code:</div><pre dir="ltr">$ bbmerge.sh in=reads.fq outa=adapters.fa reads=1m</pre></div><p><span>Of course, it will only work for paired reads! The output fasta file will look like this:</span><br /><br /></p><div><div>Code:</div><pre dir="ltr">&gt;Read1_adapter
GATCGGAAGAGCACACGTCTGAACTCCAGTCACATCACGATCTCGTATGCCGTCTTCTGCTTG
&gt;Read2_adapter
GATCGGAAGAGCACACGTCTGAACTCCAGTCACCGATGTATCTCGTATGCCGTCTTCTGCTTG</pre></div><p><span>If you have multiplexed things with different barcodes in the adapters, the part with the barcode will show up as Ns, like this:</span><br /><br /><span>GATCGGAAGAGCACACGTCTGAACTCCAGTCACNNNNNNATCTCGTATGCCGTCTTCTGCTTG&nbsp;&nbsp;</span><br /><br /><span>Note: For BBMerge with micro-RNA, you need to add the flag&nbsp;</span><strong>mininsert=17</strong><span>. The default is 35, which is too long for micro-RNA libraries.&nbsp;</span></p><ul>
<li><strong>Identifying adapters</strong></li>
</ul><p><span>If you have paired reads, and enough of the reads have inserts shorter than read length, you can identify adapter sequences with BBMerge, like this (they will be printed to adapters.fa):</span><br /><br /></p><div><div>Code:</div><pre dir="ltr">$ bbmerge.sh in1=r1.fq in2=r2.fq outa=adapters.fa</pre></div><p><br /><span>-----------------------------------------------------------------</span><br /><br /><strong>BBDuk.sh</strong><br /><br /><span>Note: BBDuk is strictly deterministic on a per-read basis, however it does by default reorder the reads when run multithreaded. You can add the flag "ordered" to keep output reads in the same order as input reads</span></p><ul>
<li><strong>Finding reads with a specific sequence at the beginning of read</strong></li>
</ul><p>&nbsp;</p><div><div>Code:</div><pre dir="ltr">$ bbduk.sh -Xmx1g in=reads.fq outm=matched.fq outu=unmatched.fq restrictleft=25 k=25 literal=AAAAACCCCCTTTTTGGGGGAAAAA</pre></div><p><span>In this case, all reads starting with "AAAAACCCCCTTTTTGGGGGAAAAA" will end up in "matched.fq" and all other reads will end up in "unmatched.fq". Specifically, the command means "look for 25-mers in the leftmost 25 bp of the read", which will require an exact prefix match, though you can relax that if you want.</span><br /><br /><span>So you could bin all the reads with your known sequence, then look at the remaining reads to see what they have in common. You can do the same thing with the tail of the read using "restrictright" instead, though you can't use both restrictions at the same time.&nbsp;&nbsp;</span><br /><br /></p><div><div>Code:</div><pre dir="ltr">$ bbduk.sh in=reads.fq outm=matched.fq literal=NNNNNNCCCCGGGGGTTTTTAAAAA k=25 copyundefined</pre></div><p><span>With the "copyundefined" flag, a copy of each reference sequence will be made representing every valid combination of defined letter. So instead of increasing memory or time use by 6^75, it only increases them by 4^6 or 4096 which is completely reasonable, but it only allows substitutions at predefined locations. You can use the "copyundefined", "hdist", and "qhdist" flags together for a lot of flexibility - for example, hdist=2 qhdist=1 and 3 Ns in the reference would allow a hamming distance of 6 with much lower resource requirements than hdist=6. Just be sure to give BBDuk as much memory as possible.</span></p><ul>
<li><strong>Removing illumina adapters (if exact adapters not known)</strong></li>
</ul><p><br /><span>If you're not sure which adapters are used, you can add "ref=truseq.fa.gz,truseq_rna.fa.gz,nextera.fa.gz" and get them all (this will increase the amount of overtrimming, though it should still be negligible).&nbsp;</span></p><ul>
<li><strong>Removing illumina control sequences/phiX reads</strong></li>
</ul><p>&nbsp;</p><div><div>Code:</div><pre dir="ltr">bbduk.sh in=trimmed.fq.gz out=filtered.fq.gz k=31 ref=artifacts,phix ordered cardinality</pre></div><ul>
<li><strong>Identify certain reads that contain a specific sequence</strong></li>
</ul><div><div>Code:</div><pre dir="ltr">$ bbduk.sh in=reads.fq out=unmatched.fq outm=matched.fq literal=ACGTACGTACGTACGTAC k=18 mm=f hdist=2</pre></div><p><span>Make sure "k" is set to the exact length of the sequence. "hdist" controls the number of substitutions allowed. "outm" gets the reads that match. By default this also looks for the reverse-complement; you can disable that with "rcomp=f".&nbsp;&nbsp;</span></p><ul>
<li><strong>Extract sequences that share kmers with your sequences with BBDuk</strong></li>
</ul><div><div>Code:</div><pre dir="ltr">$ bbduk.sh in=a.fa ref=b.fa out=c.fa mkf=1 mm=f k=31</pre></div><p><span>This will print to C all the sequences in A that share 100% of their 31-mers with sequences in B.&nbsp;</span><br /><br /></p><ul>
<li><strong>Extract sequences that contain N's with BBDuk</strong></li>
</ul><div><div>Code:</div><pre dir="ltr">bbduk.sh in=reads.fq out=readsWithoutNs.fq outm=readsWithNs.fq maxns=0</pre></div><p><span>If you have, say, 100bp reads and only want to separate reads containing all 100 Ns, change that to "maxns=99".</span><br /><br /><strong>General notes for BBDuk.sh</strong><span>&nbsp;</span><br /><br /><span>BBDuk can operate in one of 4 kmer-matching modes:</span><br /><span>Right-trimming (ktrim=r), left-trimming (ktrim=l), masking (ktrim=n), and filtering (default). But it can only do one at a time because all kmers are stored in a single table. It can still do non-kmer-based operations such as quality trimming at the same time.</span><br /><br /><span>BBDuk2 can do all 4 kmer operations at once and is designed for integration into automated pipelines where you do contaminant removal and adapter-trimming in a single pass to minimize filesystem I/O. Personally, I never use BBDuk2 from the command line. Both have identical capabilities and functionality otherwise, but the syntax is different.</span><br /><br /><span>------------------------------------------------------------------</span><br /><br /><strong>Randomreads.sh</strong></p><ul>
<li><strong>Generate random reads in various formats</strong></li>
</ul><p>&nbsp;</p><div><div>Code:</div><pre dir="ltr">$ randomreads.sh ref=genome.fasta out=reads.fq len=100 reads=10000</pre></div><p><span>You can specify paired reads, an insert size distribution, read lengths (or length ranges), and so forth. But because I developed it to benchmark mapping algorithms, it is specifically designed to give excellent control over mutations. You can specify the number of snps, insertions, deletions, and Ns per read, either exactly or probabilistically; the lengths of these events is individually customizable, the quality values can alternately be set to allow errors to be generated on the basis of quality; there's a PacBio error model; and all of the reads are annotated with their genomic origin, so you will know the correct answer when mapping.</span><br /><br /><span>Bear in mind that 50% of the reads are going to be generated from the plus strand and 50% from the minus strand. So, either a read will match the reference perfectly, OR its reverse-complement will match perfectly.</span><br /><br /><span>You can generate the same set of reads with and without SNPs by fixing the seed to a positive number, like this:</span><br /><br /></p><div><div>Code:</div><pre dir="ltr">$ randomreads.sh maxsnps=0 adderrors=false out=perfect.fastq reads=1000 minlength=18 maxlength=55 seed=5

$ randomreads.sh maxsnps=2 snprate=1 adderrors=false out=2snps.fastq reads=1000 minlength=18 maxlength=55 seed=5</pre></div><p><span>[As of BBmap v. 36.59] rendomreads.sh gains the ability to simulate metagenomes.&nbsp;</span><br /><br /><span>coverage=X will automatically set "reads" to a level that will give X average coverage (decimal point is allowed).</span><br /><br /><span>metagenome will assign each scaffold a random exponential variable, which decides the probability that a read be generated from that scaffold. So, if you concatenate together 20 bacterial genomes, you can run randomreads and get a metagenomic-like distribution. It could also be used for RNA-seq when using a transcriptome reference.</span><br /><br /><span>The coverage is decided on a per-reference-sequence level, so if a bacterial assembly has more than one contig, you may want to glue them together first with fuse.sh before concatenating them with the other references.&nbsp;</span><br /><br /></p><ul>
<li><strong>Simulate a jump library</strong></li>
</ul><p><br /><span>You can simulate a 4000bp jump library from your existing data like this.</span><br /><br /></p><div><div>Code:</div><pre dir="ltr">$ cat assembly1.fa assembly2.fa &gt; combined.fa
$ bbmap.sh ref=combined.fa
$ randomreads.sh reads=1000000 length=100 paired interleaved mininsert=3500 maxinsert=4500 bell perfect=1 q=35 out=jump.fq.gz</pre></div><p><span>--------------------------------------------------------------</span><br /><strong>Shred.sh</strong><br /><br /></p><div><div>Code:</div><pre dir="ltr">$ shred.sh in=ref.fasta out=reads.fastq length=200</pre></div><p><span>The difference is that RandomReads will make reads in a random order from random locations, ensuring flat coverage on average, but it won't ensure 100% coverage unless you generate many fold depth. Shred, on the other hand, gives you exactly 1x depth and exactly 100% coverage (and is not capable of modelling errors). So, the use-cases are different.&nbsp;</span><br /><span>---------------------------------------------------------------</span><br /><strong>Demuxbyname.sh</strong></p><ul>
<li><strong>Demultiplex fastq files when the tag is present in the fastq read header (illumina)</strong></li>
</ul><p>&nbsp;</p><div><div>Code:</div><pre dir="ltr">$ demuxbyname.sh in=r#.fq out=out_%_#.fq prefixmode=f names=GGACTCCT+GCGATCTA,TAAGGCGA+TCTACTCT,...
outu=filename</pre></div><p><span>"Names" can also be a text file with one barcode per line (in exactly the format found in the read header). You do have to include all of the expected barcodes, though.</span><br /><br /><span>In the output filename, the "%" symbol gets replaced by the barcode; in both the input and output names, the "#" symbol gets replaced by 1 or 2 for read 1 or read 2. It's optional, though; you can leave it out for interleaved input/output, or specify in1=/in2=/out1=/out2= if you want custom naming.</span><br /><br /><span>----------------------------------------------------------------</span><br /><br /><strong>Readlength.sh</strong></p><ul>
<li><strong>Plotting the length distribution of reads</strong></li>
</ul><div><div>Code:</div><pre dir="ltr">$ readlength.sh in=file out=histogram.txt bin=10 max=80000</pre></div><p><span>That will plot the result in bins of size 10, with everything above 80k placed in the same bin. The defaults are set for relatively short sequences so if they are many megabases long you may need to add the flag "-Xmx8g" and increase "max=" to something much higher.</span><br /><br /><span>Alternatively, if these are assemblies and you're interested in continuity information (L50, N50, etc), you can run stats on each or statswrapper on all of them:</span><br /><br /></p><div><div>Code:</div><pre dir="ltr">stats.sh in=file</pre></div><p><span>or</span><br /><br /></p><div><div>Code:</div><pre dir="ltr">statswrapper.sh in=file,file,file,file&hellip;</pre></div><p><span>----------------------------------------------------------------</span><br /><strong>Filterbyname.sh</strong><br /><br /><span>By default, "filterbyname" discards reads with names in your name list, and keeps the rest. To include them and discard the others, do this:</span><br /><br /></p><div><div>Code:</div><pre dir="ltr">$ filterbyname.sh in=003.fastq out=filter003.fq names=names003.txt include=t</pre></div><p><span>----------------------------------------------------------------</span><br /><strong>getreads.sh</strong><br /><br /><span>If you only know the number(s) of the fasta/fastq record(s) in a file (records start at 0) then you can use the following command to extract those reads in a new file.</span><br /><br /></p><div><div>Code:</div><pre dir="ltr">$ getreads.sh in= id=&lt;number,number,number...&gt; out=</pre></div><p><span>The first read (or pair) has ID 0, the second read (or pair) has ID 1, etc.</span><br /><br /><span>Parameters:</span><br /><span>in= Specify the input file, or stdin.</span><br /><span>out= Specify the output file, or stdout.</span><br /><span>id= Comma delimited list of numbers or ranges, in any order.</span><br /><span>For example: id=5,93,17-31,8,0,12-13&nbsp;</span><br /><span>----------------------------------------------------------------</span><br /><strong>Splitsam.sh</strong></p><ul>
<li><strong>Splits a sam file into forward and reverse reads</strong></li>
</ul><p>&nbsp;</p><div><div>Code:</div><pre dir="ltr">splitsam.sh mapped.sam plus.sam minus.sam unmapped.sam
reformat.sh in=plus.sam out=plus.fq
reformat.sh in=minus.sam out=minus.fq rcomp</pre></div><p><span>----------------------------------------------------------------</span><br /><strong>BBSplit.sh</strong><br /><br /><span>BBSplit now has the ability to output paired reads in dual files using the # symbol. For example:</span><br /><br /></p><div><div>Code:</div><pre dir="ltr">$ bbsplit.sh ref=x.fa,y.fa in1=read1.fq in2=read2.fq basename=o%_#.fq</pre></div><p><span>will produce ox_1.fq, ox_2.fq, oy_1.fq, and oy_2.fq</span><br /><br /><span>You can use the # symbol for input also, like "in=read#.fq", and it will get expanded into 1 and 2.&nbsp;&nbsp;</span><br /><br /><strong>Added feature:&nbsp;</strong><span>One can specify a directory for the "ref=" argument. If anything in the list is a directory, it will use all fasta files in that directory. They need a fasta extension, like .fa or .fasta, but can be compressed with an additional .gz after that. Reason this is useful is to use BBSplit is to have it split input into one output file per reference file.</span><br /><br /><br /><strong>NOTE: 1</strong><span>&nbsp;By default BBSplit uses fairly strict mapping parameters; you can get the same sensitivity as BBMap by adding the flags "minid=0.76 maxindel=16k minhits=1". With those parameters it is extremely sensitive.</span><br /><br /><strong>NOTE: 2</strong><span>&nbsp;BBSplit has different ambiguity settings for dealing with reads that map to multiple genomes. In any case, if the alignment score is higher to one genome than another, it will be associated with that genome only (this considers the combined scores of read pairs - pairs are always kept together). But when a read or pair has two identically-scoring mapping locations, on different genomes, the behavior is controlled by the "ambig2" flag - "ambig2=toss" will discard the read, "all" will send it to all output files, and "split" will send it to a separate file for ambiguously-mapped reads (one per genome to which it maps).</span><br /><br /><strong>NOTE: 3</strong><span>&nbsp;Zero-count lines are suppressed by default, but they should be printed if you include the flag "nzo=f" (nonzeroonly=false).&nbsp;</span><br /><br /><strong>NOTE: 4</strong><span>&nbsp;BBSplit needs multiple reference files as input; one per organism, or one for target and another for everything else. It only outputs one file per reference file.</span><br /><br /><span>Seal.sh, on the other hand, which is similar, can use a single concatenated file, as it (by default) will output one file per reference sequence within a concatenated set of references.&nbsp;</span><br /><span>--------------------------------------------------------------</span><br /><strong>Pileup.sh</strong></p><ul>
<li><strong>To generate transcript coverage stats</strong></li>
</ul><p>&nbsp;</p><div><div>Code:</div><pre dir="ltr">$ pileup.sh in=mapped.sam normcov=normcoverage.txt normb=20 stats=stats.txt</pre></div><p><span>That will generate coverage per transcript, with 20 lines per transcript, each line showing the coverage for that fraction of the transcript. "stats" will contain other information like the fraction of bases in each transcript that was covered.&nbsp;</span></p><ul>
<li><strong>To calculate physical coverage stats (region covered by paired-end reads)&nbsp;</strong></li>
</ul><p><span>BBMap has a "physcov" flag that allows it to report physical rather than sequenced coverage. It can be used directly in BBMap, or with pileup, if you already have a sam file. For example:</span><br /><br /></p><div><div>Code:</div><pre dir="ltr">$ pileup.sh in=mapped.sam covstats=coverage.txt</pre></div><ul>
<li><strong>Calculating coverage of the genome</strong></li>
</ul><p><br /><span>Program will take sam or bam, sorted or unsorted.</span><br /><br /></p><div><div>Code:</div><pre dir="ltr">$ pileup.sh in=mapped.sam out=stats.txt hist=histogram.txt</pre></div><p><span>stats.txt will contain the average depth and percent covered of each reference sequence; the histogram will contain the exact number of bases with a each coverage level. You can also get per-base coverage or binned coverage if you want to plot the coverage. It also generates median and standard deviation, and so forth.</span><br /><br /><span>It's also possible to generate coverage directly from BBMap, without an intermediate sam file, like this:</span><br /><br /></p><div><div>Code:</div><pre dir="ltr">$ bbmap.sh in=reads.fq ref=reference.fasta nodisk covstats=stats.txt covhist=histogram.txt</pre></div><p><span>We use this a lot in situations where all you care about is coverage distributions, which is somewhat common in metagenome assemblies. It also supports most of the flags that pileup.sh supports, though the syntax is slightly different to prevent collisions. In each case you can see all the possible flags by running the shellscript with no arguments.</span></p><ul>
<li><strong>To bin aligned reads</strong></li>
</ul><p>&nbsp;</p><div><div>Code:</div><pre dir="ltr">$ pileup.sh in=mapped.sam out=stats.txt bincov=coverage.txt binsize=1000</pre></div><p><span>That will give coverage within each bin. For read density regardless of read length, add the "startcov=t" flag.&nbsp;&nbsp;</span><br /><br /><span>--------------------------------------------------------------</span><br /><strong>Dedupe.sh</strong><br /><br /><span>Dedupe ensures that there is at most one copy of any input sequence, optionally allowing contaminants (substrings) to be removed, and a variable hamming or edit distance to be specified. Usage:</span><br /><br /></p><div><div>Code:</div><pre dir="ltr">$ dedupe.sh in=assembly1.fa,assembly2.fa out=merged.fa</pre></div><p><span>That will absorb exact duplicates and containments. You can use "hdist" and "edist" flags to allow mismatches, or get a complete list of flags by running the shellscript with no arguments.&nbsp;&nbsp;</span><br /><br /><span>Dedupe&nbsp;</span><span style="text-decoration: underline;">will merge assemblies</span><span>, but it&nbsp;</span><span style="text-decoration: underline;">will not produce consensus sequences or join overlapping reads</span><span>; it only removes sequences that are fully contained within other sequences (allowing the specified number of mismatches or edits).</span><br /><br /><span>Dedupe can remove duplicate reads from multiple files simultaneously, if they are comma-delimited (e.g. in=file1.fastq,file2.fastq,file3.fastq). And if you set the flag "uniqueonly=t" then ALL copies of duplicate reads will be removed, as opposed to the default behavior of leaving one copy of duplicate reads.</span><br /><br /><span>However, it does not care which file a read came from; in other words, it can't remove only reads that are duplicates across multiple files but leave the ones that are duplicates within a file. That can still be accomplished, though, like this:</span><br /><br /><span>1) Run dedupe on each sample individually, so now there are at most 1 copy of a read per sample.</span><br /><span>2) Run dedupe again on all of the samples together, with "uniqueonly=t". The only remaining duplicate reads will be the ones duplicated between samples, so that's all that will be removed.&nbsp;&nbsp;</span><br /><br /><span>--------------------------------------------------------------</span></p><ul>
<li><strong>Generate ROC curves from any aligner</strong></li>
</ul><p><br /><strong>[*]index the reference<br /><br /></strong></p><div><div>Code:</div><pre dir="ltr">$ bbmap.sh ref=reference.fasta</pre></div><p><br /><strong>[*]Generate random reads</strong><br /><br /></p><div><div>Code:</div><pre dir="ltr">$ randomreads.sh reads=100000 length=100 out=synth.fastq maxq=35 midq=25 minq=15</pre></div><p><strong>[*]Map to produce a sam file</strong><br /><br /><span>...substitute this command with the appropriate one from your aligner of choice</span></p><div><div>Code:</div><pre dir="ltr">$ bbmap.sh in=synth.fq out=mapped.sam</pre></div><p><strong>[*]Generate ROC curve</strong><br /><br /></p><div><div>Code:</div><pre dir="ltr">$ samtoroc.sh in=mapped.sam reads=100000</pre></div><p><span>--------------------------------------------------------------</span></p><ul>
<li><strong>Calculate heterozygous rate for sequence data</strong></li>
</ul><p><br /><br /></p><div><div>Code:</div><pre dir="ltr">$ kmercountexact.sh in=reads.fq khist=histogram.txt peaks=peaks.txt</pre></div><p><span>You can examine the histogram manually, or use the "peaks" file which tells you the number of unique kmers in each peak on the histogram. For a diploid, the first peak will be the het peak, the second will be the homozygous peak, and the rest will be repeat peaks. The peak caller is not perfect, though, so particularly with noisy data I would only rely on it for the first two peaks, and try to quantify the higher-order peaks manually if you need to (which you generally don't).</span><br /><br /><span>-----------------------------------------------------------------</span></p><ul>
<li><strong>Compare mapped reads between two files</strong></li>
</ul><p><br /><span>To see how many mapped reads (can be mapped concordant or discordant, doesn't matter) are shared between the two alignment files and how many mapped reads are unique to one file or the other.</span><br /><br /></p><div><div>Code:</div><pre dir="ltr">$ reformat.sh in=file1.sam out=mapped1.sam mappedonly
$ reformat.sh in=file2.sam out=mapped2.sam mappedonly</pre></div><p><span>That gets you the mapped reads only. Then:</span><br /><br /></p><div><div>Code:</div><pre dir="ltr">$ filterbyname.sh in=mapped1.sam names=mapped2.sam out=shared.sam include=t</pre></div><p><span>...which gets you the set intersection;</span><br /><br /></p><div><div>Code:</div><pre dir="ltr">$ filterbyname.sh in=mapped1.sam names=mapped2.sam out=only1.sam include=f
$ filterbyname.sh in=mapped2.sam names=mapped1.sam out=only2.sam include=f</pre></div><p><span>...which get you the set subtractions.&nbsp;&nbsp;</span><br /><br /><span>--------------------------------------------------------------</span><br /><br /><strong>BBrename.sh</strong></p><div><div>Code:</div><pre dir="ltr">$ bbrename.sh in=old.fasta out=new.fasta</pre></div><p><span>That will rename the reads as 1, 2, 3, 4, ... 222.</span><br /><br /><span>You can also give a custom prefix if you want. The input has to be text format, not .doc.&nbsp;&nbsp;</span><br /><br /><span>---------------------------------------------------------------------</span><br /><br /><strong>BBfakereads.sh</strong></p><ul>
<li><strong>Generating &ldquo;fake&rdquo; paired end reads from a single end read file</strong></li>
</ul><p>&nbsp;</p><div><div>Code:</div><pre dir="ltr">$ bfakereads.sh in=reads.fastq out1=r1.fastq out2=r2.fastq length=100</pre></div><p><span>That will generate fake pairs from the input file, with whatever length you want (maximum of input read length). We use it in some cases for generating a fake LMP library for scaffolding from a set of contigs. Read 1 will be from the left end, and read 2 will be reverse-complemented and from the right end; both will retain the correct original qualities. And " /1" " /2" will be suffixed after the read name.&nbsp;&nbsp;</span><br /><br /><span>------------------------------------------------------------------</span><br /><strong>Randomreads.sh</strong></p><ul>
<li><strong>Generate random reads</strong></li>
</ul><p>&nbsp;</p><div><div>Code:</div><pre dir="ltr">$ randomreads.sh ref=genome.fasta out=reads.fq len=100 reads=10000</pre></div><p><span>"seed=-1" will use a random seed; any other value will use that specific number as the seed</span><br /><br /><span>You can specify paired reads, an insert size distribution, read lengths (or length ranges), and so forth. But because I developed it to benchmark mapping algorithms, it is specifically designed to give excellent control over mutations. You can specify the number of snps, insertions, deletions, and Ns per read, either exactly or probabilistically; the lengths of these events is individually customizable, the quality values can alternately be set to allow errors to be generated on the basis of quality; there's a PacBio error model; and all of the reads are annotated with their genomic origin, so you will know the correct answer when mapping.</span><br /><br /><span>--------------------------------------------------------------------</span></p><ul>
<li><strong>Generate saturation curves to assess sequencing depth</strong></li>
</ul><p>&nbsp;</p><div><div>Code:</div><pre dir="ltr">$ bbcountunique.sh in=reads.fq out=histogram.txt</pre></div><p><span>It works by pulling kmers from each input read, and testing whether it has been seen before, then storing it in a table.</span><br /><br /><span>The bottom line, "first", tracks whether the first kmer of the read has been seen before (independent of whether it is read 1 or read 2).</span><br /><br /><span>The top line, "pair", indicates whether a combined kmer from both read 1 and read 2 has been seen before. The other lines are generally safe to ignore but they track other things, like read1- or read2-specific data, and random kmers versus the first kmer.</span><br /><br /><span>It plots a point every X reads (configurable, default 25000).</span><br /><br /><span>In noncumulative mode (default), a point indicates "for the last X reads, this percentage had never been seen before". In this mode, once the line hits zero, sequencing more is not useful.</span><br /><br /><span>In cumulative mode, a point indicates "for all reads, this percentage had never been seen before", but still only one point is plotted per X reads.</span><br /><br /><span>-----------------------------------------------------------------</span><br /><strong>CalcTrueQuality.sh</strong><br /><br /><a href="http://seqanswers.com/forums/showthread.php?p=170904" target="_blank">http://seqanswers.com/forums/showthread.php?p=170904</a><br /><br /><span>In light of the quality-score issues with the NextSeq platform, and the possibility of future Illumina platforms (HiSeq 3000 and 4000) also using quantized quality scores, I developed it for recalibrating the scores to ensure accuracy and restore the full range of values.</span><br /><br /><span>-----------------------------------------------------------------</span><br /><br /><strong>BBMapskimmer.sh</strong><br /><br /><span>BBMap is designed to find the best mapping, and heuristics will cause it to ignore mappings that are valid but substantially worse. Therefore, I made a different version of it, BBMapSkimmer, which is designed to find all of the mappings above a certain threshold. The shellscript is bbmapskimmer.sh and the usage is similar to bbmap.sh or mapPacBio.sh. For primers, which I assume will be short, you may wish to use a lower than default K of, say, 10 or 11, and add the "slow" flag.</span><br /><br /><span>--------------------------------------------------------------</span><br /><br /><strong>msa.sh and curprimers.sh</strong><br /><br /><span>Quoted from Brian's response directly.</span><br /><br /><span>I also wrote another pair of programs specifically for working with primer pairs, msa.sh and cutprimers.sh. msa.sh will forcibly align a primer sequence (or a set of primer sequences) against a set of reference sequences to find the single best matching location per reference sequence - in other words, if you have 3 primers and 100 ref sequences, it will output a sam file with exactly 100 alignments - one per ref sequence, using the primer sequence that matched best. Of course you can also just run it with 1 primer sequence.</span><br /><br /><span>So you run msa twice - once for the left primer, and once for the right primer - and generate 2 sam files. Then you feed those into cutprimers.sh, which will create a new fasta file containing the sequence between the primers, for each reference sequence. We used these programs to synthetically cut V4 out of full-length 16S sequences.</span><br /><br /><span>I should say, though, that the primer sites identified are based on the normal BBMap scoring, which is not necessarily the same as where the primers would bind naturally, though with highly conserved regions there should be no difference.</span><br /><br /><span>------------------------------------------------------</span><br /><strong>testformat.sh</strong><br /><br /><strong>Identify type of Q-score encoding in sequence files</strong><br /><br /></p><div><div>Code:</div><pre dir="ltr">$ testformat.sh in=seq.fq.gz
sanger    fastq    gz    interleaved    150bp</pre></div><p><span>--------------------------------------------------</span><br /><strong>kcompress.sh</strong><br /><br /><span>Newest member of BBTools. Identify constituent k-mers.&nbsp;</span><br /><a href="http://seqanswers.com/forums/showthread.php?t=63258" target="_blank">http://seqanswers.com/forums/showthread.php?t=63258</a><br /><br /><span>----------------------------------------------------</span><br /><strong>commonkmers.sh</strong><br /><br /><span>Find all k-mers for a given sequence.</span></p><div><div>Code:</div><pre dir="ltr">$ commonkmers.sh in=reads.fq out=kmers.txt k=4 count=t display=999</pre></div><p><span>Will produce output that looks like</span><br /><br /></p><div><div>Code:</div><pre dir="ltr">MISEQ05:239:000000000-A74HF:1:2110:14788:23085	ATGA=8	ATGC=6	GTCA=6	AAAT=5	AAGC=5	AATG=5	AGCA=5	ATAA=5	ATTA=5	CAAA=5	CATA=5	CATC=5	CTGC=5	AACC=4	AACG=4	AAGA=4	ACAT=4	ACCA=4	AGAA=4	ATCA=4	ATGG=4	CAAG=4	CCAA=4	CCTC=4	CTCA=4	CTGA=4	CTTC=4	GAGC=4	GGTA=4	GTAA=4	GTTA=4	AAAA=3	AAAC=3	AAGT=3	ACCG=3	ACGG=3	ACTG=3	AGAT=3	AGCT=3	AGGA=3	AGTA=3	AGTC=3	CAGC=3	CATG=3	CGAG=3	CGGA=3	CGTC=3	CTAA=3	CTCC=3	CTTA=3	GAAA=3	GACA=3	GACC=3	GAGA=3	GCAA=3	GGAC=3	TCAA=3	TGCA=3	AAAG=2	AACA=2	AATA=2	AATC=2	ACAA=2	ACCC=2	ACCT=2	ACGA=2	ACGC=2	AGAC=2	AGCG=2	AGGC=2	CAAC=2	CAGG=2	CCGC=2	GCCA=2	GCTA=2	GGAA=2	GGCA=2	TAAA=2	TAGA=2	TCCA=2	TGAA=2	AAGG=1	AATT=1	ACGT=1	AGAG=1	AGCC=1	AGGG=1	ATAC=1	ATAG=1	ATTG=1	CACA=1	CACG=1	CAGA=1	CCAC=1	CCCA=1	CCGA=1	CCTA=1	CGAC=1	CGCA=1	CGCC=1	CGCG=1	CGTA=1	CTAC=1	GAAC=1	GCGA=1	GCGC=1	GTAC=1	GTGA=1	TTAA=1</pre></div><p><span>-----------------------------------------------------</span><br /><strong>Mutate.sh</strong><br /><br /><span>Simulate multiple mutants from a known reference (e.g.&nbsp;</span><em>E. coli</em><span>).</span><br /><br /></p><div><div>Code:</div><pre dir="ltr">$ mutate.sh in=e_coli.fasta out=mutant.fasta id=99 
$ randomreads.sh ref=mutant.fasta out=reads.fq.gz reads=5m length=150 paired adderrors</pre></div><p><span>That will create a mutant version of E.coli with 99% identity to the original, and then generate 5 million simulated read pairs from the new genome. You can repeat this multiple times; each mutant will be different.</span><br /><br /><span>------------------------------------</span><br /><br /><strong>Partition.sh</strong><br /><br /><span>One can partition a large dataset with partition.sh into smaller subsets (example below splits data into 8 chunks).</span><br /><br /></p><div><div>Code:</div><pre dir="ltr">partition.sh in=r1.fq in2=r2.fq out=r1_part%.fq out2=r2_part%.fq ways=8</pre></div><p><span>-----------------------------------</span><br /><strong>clumpify.sh</strong><br /><br /><span>If you are concerned about file size and want the files to be as small as possible, give Clumpify a try. It can reduce filesize by around 30% losslessly by reordering the reads. I've found that this also typically accelerates subsequent analysis pipelines by a similar factor (up to 30%). Usage:</span><br /><br /></p><div><div>Code:</div><pre dir="ltr">clumpify.sh in=reads.fastq.gz out=clumped.fastq.gz</pre></div><div><div>Code:</div><pre dir="ltr">clumpify.sh in1=reads_R1.fastq.gz in2=reads_R2.fastq.gz out1=clumped_R1.fastq.gz out2=clumped_R2.fastq.gz</pre></div><ul>
<li><strong>Clumpify.sh can now mark/remove sequence duplicates (optical/PCR/otherwise) from NGS data</strong></li>
</ul><p><br /><span>This does NOT require alignments so it should prove more useful compared to Picard MarkDuplicates. Relevant options for clumpify.sh command are listed below.</span><br /><br /></p><div><div>Code:</div><pre dir="ltr">dedupe=f optical=f (default)
Nothing happens with regards to duplicates.

dedupe=t optical=f
All duplicates are detected, whether optical or not.  All copies except one are removed for each duplicate.

dedupe=f optical=t
Nothing happens.

dedupe=t optical=t

Only optical duplicates (those with an X or Y coordinate within dist) are detected.  All copies except one are removed for each duplicate.
The allduplicates flag makes all copies of duplicates removed, rather than leaving a single copy.  But like optical, it has no effect unless dedupe=t.

Note: If you set "dupedist" to anything greater than 0, "optical" gets enabled automatically.</pre></div><p><span>-------------------------------------</span><br /><strong>fuse.sh</strong><br /><br /><span>Fuse will automatically reverse-complement read 2. Pad (N) amount can be adjusted as necessary. This will for example create a full size amplicon that can be used for alignments.</span><br /><br /></p><div><div>Code:</div><pre dir="ltr">fuse.sh in1=r1.fq in2=r2.fq pad=130 out=fused.fq fusepairs</pre></div>]]></description>
	<dc:creator>Surabhi Chaudhary</dc:creator>
</item>
<item>
	<guid isPermaLink="true">https://bioinformaticsonline.com/blog/view/35376/how-custom-essay-writing-service-work-to-provide-essay-papers-that-meets-customers-need</guid>
	<pubDate>Sat, 27 Jan 2018 04:13:32 -0600</pubDate>
	<link>https://bioinformaticsonline.com/blog/view/35376/how-custom-essay-writing-service-work-to-provide-essay-papers-that-meets-customers-need</link>
	<title><![CDATA[How custom essay writing service work to provide essay papers that meets customers need?]]></title>
	<description><![CDATA[<p>Writing is more effective way to learn than reading. While writing you get more information and hence can improve your personal skills and knowledge. Every student likes to get a good essay, for getting good essay students must follow the some procedures to write essay. There are many more online writing services available to help students who can&rsquo;t complete their writing papers.</p><p>&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp; Custom essay writing service provides better service than other service available. Because, here the working procedure is very simple. It is a very easy task for the customer to get an essay with this site because the process is straightforward. There are just four steps that the customer needs to follow. The first is to place the order. There are many writers who are specialized to their subject. This writer will carefully read through the requirements and give a rough draft of the final paper and will be one or two pages long. Which is the example ,that the final essay paper look like. This paper is then send to the customer by the writers. The communication between the writer and the customer is very important to pass information to the writer any time of the day. The final step is the submission of the paper ,also provide time for the paper to be read by the customer and any revisions will have to be made before submission. The customer can do any type of revision request and the writer will do any number of revisions till the customer is satisfied with the paper.</p><p>&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp; Hence we can conclude that <a href="https://essaycaptains.com/">custom essay writing service</a> work for customer at any time of the day and give a better result to the customers.</p><p>&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;</p>]]></description>
	<dc:creator>melissamfloyd</dc:creator>
</item>
<item>
	<guid isPermaLink="true">https://bioinformaticsonline.com/blog/view/34916/bioinformatics-tools-developed-for-oxford-nanopore-data-analysis</guid>
	<pubDate>Wed, 27 Dec 2017 20:47:30 -0600</pubDate>
	<link>https://bioinformaticsonline.com/blog/view/34916/bioinformatics-tools-developed-for-oxford-nanopore-data-analysis</link>
	<title><![CDATA[Bioinformatics tools developed for Oxford Nanopore data analysis !]]></title>
	<description><![CDATA[<p><span>MinION is the only portable real-time device for DNA and RNA&nbsp;</span><span>sequencing</span><span>. Each consumable flow cell can now generate 10&ndash;20 Gb of DNA&nbsp;</span><span>sequence</span><span>&nbsp;data. Ultra-</span><span>long read lengths are possible (hundreds of kb) as you can choose your fragment length.&nbsp;</span>One of the technical advantages of ONT data is the read length, which offers great prospects for genome assembly. Generally, assemblers are based on several different types of algorithms, such as greedy, overlap-layout-consensus (OLC), de Bruijn graph (DBG), and string graph.</p><p><span>List of analysis tools developed for Oxford Nanopore data</span></p><p>BWA <br />Fast nanopore data tuned alignment tool <br />https://github.com/lh3/bwa</p><p>GraphMap<br />Mapper for long and error-prone reads<br />https://github.com/isovic/graphmap</p><p>LAST<br />Nanopore tuned alignment tool<br />http://last.cbrc.jp/</p><p>LINKS<br />Software tool for long read scaffolding <br />https://github.com/warrenlr/LINKS/</p><p>marginAlign<br />Tools to align nanopore reads to a reference<br />https://github.com/benedictpaten/marginAlign</p><p>minoTour<br />Real time analysis tools<br />http://minotour.nottingham.ac.uk/</p><p>nanoCORR<br />Error-correction tool for nanopore sequence data<br />https://github.com/jgurtowski/nanocorr</p><p>NanoOK<br />Software for nanopore data, quality and error profiles<br />https://documentation.tgac.ac.uk/display/NANOOK/NanoOK</p><p>Nanopolish<br />Nanopore analysis and genome assembly software<br />https://github.com/jts/nanopolish</p><p>nanopore<br />Variant-detection tool for nanopore sequence data<br />https://github.com/mitenjain/nanopore</p><p>Nanocorrect<br />Error-correction tool for nanopore sequence data<br />https://github.com/jts/nanocorrect/</p><p>npReader<br />Real-time conversion and analysis of nanopore reads<br />https://github.com/mdcao/npReader</p><p>poRe<br />Tool for analyzing and visualizing nanopore data<br />https://sourceforge.net/p/rpore/wiki/Home/</p><p>PoreSeq<br />Error-correction and variant-calling software<br />https://github.com/tszalay/poreseq</p><p>Poretools<br />Nanopore sequence analysis and visualization software <br />https://github.com/arq5x/poretools</p><p>SSPACE-LongRead<br />Genome scaffolding tool <br />http://www.baseclear.com/genomics/bioinformatics/basetools/SSPACE-longread</p><p>SMIS<br />Genome scaffolding tool <br />https://sourceforge.net/projects/phusion2/files/smis/</p><p>&nbsp;</p><p>List of assemblers for Oxford Nanopore MinION long reads</p><p>LQS<br />DALIGNER, Celera OLC Nanocorrect, <br />Nanopolish corrector<br />https://github.com/jts/nanopolish</p><p>PBcR<br />HGAP or BLASR, Celera OLC <br />PBcR corrector<br />http://wgs-assembler.sourceforge.net/wiki/index.php/PBcR<br /> &ndash;<br />Canu<br />MHAP, Celera OLC <br />Canu corrector<br />https://github.com/marbl/canu</p><p>Falcon<br />String graph, Celera OLC <br />Falcon corrector<br />https://github.com/PacificBiosciences/falcon</p><p>Miniasm <br />OLC<br />https://github.com/lh3/miniasm</p><p>ra-integrate<br />OLC<br />https://github.com/mariokostelac/ra-integrate/</p><p>ALLPATHS-LG<br />de Bruijn graph <br />ALLPATHS-L corrector<br />https://www.broadinstitute.org/software/allpaths-lg/blog/?page_id=12</p><p>SPAdes <br />de Bruijn graph <br />SPAdes corrector<br />http://bioinf.spbau.ru/spades</p>]]></description>
	<dc:creator>biogeek</dc:creator>
</item>
<item>
	<guid isPermaLink="true">https://bioinformaticsonline.com/blog/view/34912/list-of-cancer-genomics-research-web-resources</guid>
	<pubDate>Wed, 27 Dec 2017 20:33:09 -0600</pubDate>
	<link>https://bioinformaticsonline.com/blog/view/34912/list-of-cancer-genomics-research-web-resources</link>
	<title><![CDATA[List of cancer genomics research web resources !]]></title>
	<description><![CDATA[<p>Major web resources for cancer genomics research</p><p>CGHub <br />https://cghub.ucsc.edu/ <br />Comprehensive data repository; huge data size</p><p>EGA <br />https://www.ebi.ac.uk/ega/ <br />Comprehensive data repository; huge data size</p><p>COSMIC <br />http://cancer.sanger.ac.uk <br />Largest somatic mutation database; genome sequencing paper curation</p><p>CPRG <br />http://www.broadinstitute.org/software/cprg <br />Interface for cancer program resources</p><p>GDAC <br />http://gdac.broadinstitute.org/ <br />Data analysis; automatic pipelines; user-friendly reports</p><p>SNP500Cancer <br />http://snp500cancer.nci.nih.gov <br />Sequence and genotype verification of SNPs</p><p>canEvolve <br />www.canevolve.org/ <br />Comprehensive analysis of tumor profile; Data from 90 studies involving more than 10,000 patients</p><p>MethyCancer <br />http://methycancer.psych.ac.cn <br />Relationship among DNA methylation, gene expression and cancer</p><p>SomamiR <br />http://compbio.uthsc.edu/SomamiR/ <br />Correlation between somatic mutation and microRNA; genome-wide displaying</p><p>cBioPortal <br />http://www.cbioportal.org/public-portal/ <br />Graphical summaries; gene alteration; processed data; visualization</p><p>UCSC Cancer Genomics Browser <br />https://genome-cancer.soe.ucsc.edu/ <br />Clinical information; gene expression; copy number variation; visualization</p><p>CGWB <br />https://cgwb.nci.nih.gov/ <br />Visualization; gene mutation and variation; automated analysis pipeline</p><p>GDSC <br />http://www.cancerrxgene.org <br />Drug sensitivity information; drug response information</p><p>canSAR <br />https://cansar.icr.ac.uk/ <br />Multidisciplinary information; drug discovery</p><p>NONCODE <br />http://www.noncode.org/ ncRNAs; <br />lncRNAs; up-to-date and comprehensive resource</p>]]></description>
	<dc:creator>biogeek</dc:creator>
</item>
<item>
	<guid isPermaLink="true">https://bioinformaticsonline.com/blog/view/34864/installing-perl-environment-on-linux</guid>
	<pubDate>Tue, 26 Dec 2017 21:21:50 -0600</pubDate>
	<link>https://bioinformaticsonline.com/blog/view/34864/installing-perl-environment-on-linux</link>
	<title><![CDATA[Installing Perl environment on Linux]]></title>
	<description><![CDATA[<p>By using&nbsp;<code>plenv</code>, you can easily install and switch among different version of Perl. This will be installed under your home directory in<code>~/.plenv</code>.</p><h4>Install latest Perl (with supporting multithreading) and CPANMinus.</h4><pre><code> $ cd
 $ git clone git://github.com/tokuhirom/plenv.git ~/.plenv
 $ git clone git://github.com/tokuhirom/Perl-Build.git ~/.plenv/plugins/perl-build/
 $ echo 'export PATH="$HOME/.plenv/bin:$PATH"' &gt;&gt; ~/.bashrc
 $ echo 'eval "$(plenv init -)"' &gt;&gt; ~/.bashrc
 $ source ~/.bashrc
 $ plenv install 5.18.1 -Dusethreads
 $ plenv rehash
 $ plenv global 5.18.1
 $ plenv install-cpanm
</code></pre><ul>
<li><code>git</code>&nbsp;is a distributed revision control and source code management software which can help you to download files from GitHub server.</li>
<li><code>echo</code>&nbsp;means "print".</li>
<li><code>&gt;&gt;</code>&nbsp;means adding the output into the end of the file, while&nbsp;<code>&gt;</code>&nbsp;means adding the output by overwriting the whole file. Please use<code>&gt;</code>&nbsp;with additional cares.</li>
<li>In Linux system, there are two types of outputs when you execute a command. One is called standard output (or sometimes STDOUT for short), and the other is a standard error (STDERR).&nbsp;<code>1&gt;</code>&nbsp;is for STDOUT only,&nbsp;<code>2&gt;</code>&nbsp;is for STDERR only, and&nbsp;<code>&amp;&gt;</code>means for both. In default&nbsp;<code>&gt;</code>&nbsp;is the same to&nbsp;<code>1&gt;</code>.</li>
<li><code>exec</code>&nbsp;is execution.</li>
<li>Remember to install Perl in supporting multithreading (with option&nbsp;<code>-Dusethreads</code>), which is important for many NGS analysis packages (e.g. Trinity). In this setting, you can use multiple CPU for Perl software.</li>
<li>Install the CPAN (Comprehensive Perl Archive Network) manager software, CPANMinus, by&nbsp;<code>install-cpanm</code>.</li>
</ul><p>You can use&nbsp;<code>plenv global</code>&nbsp;and&nbsp;<code>plenv local</code>&nbsp;to change the different version of Perl to fulfil different needs of your Perl software.</p><p>For example, if the&nbsp;specific version of Perl is not compatible with your script, you can switch to the different version by:</p><pre><code> $ plenv local 
</code></pre><ul>
<li>It is similar to set the local version of your script language when you use&nbsp;<code>pyenv</code>&nbsp;and&nbsp;<code>rbenv</code>&nbsp;as the following.</li>
</ul><p>Put the following path into&nbsp;<code>~/.bashrc file</code>.</p><pre><code>export PERL5LIB="$HOME/.plenv/build/perl-5.18.1/lib"
</code></pre><h4>Install BioPerl and PerlIO::gzip</h4><p>CPANMinus is a very good Perl module manager, use&nbsp;<code>cpanm</code>&nbsp;to install BioPerl can save you a lot of time. Here are some useful modules:</p><pre><code>$ cpanm Bio::Perl
$ cpanm Bio::SearchIO
$ cpanm PerlIO::gzip<br /></code></pre><p><span>For more information, please visit:&nbsp;</span><a href="https://github.com/tokuhirom/plenv">https://github.com/tokuhirom/plenv</a></p><pre><code>&nbsp;</code></pre>]]></description>
	<dc:creator>biogeek</dc:creator>
</item>

</channel>
</rss>