Internals
---------

Overview
~~~~~~~~

.. figure:: http://chapmanb.github.com/bcbb/nglims_organization.png
   :align: center
   :alt: Overview of workflow

   Overview of workflow
   
Code structure
~~~~~~~~~~~~~~

The main scripts that handle automation of the analysis and storage are:

-  ``scripts/illumina_finished_msg.py`` -- Sits on the sequencer output
   machine; run via a cron job every hour to check for new output runs.

-  ``scripts/nextgen_analysis_server.py`` -- Main server script; runs
   specific servers for top level analysis management, storage, or
   distributed processing of individual analysis steps:

-  Called with ``-q toplevel`` -- Coordinate sample processing, running
   the full automated analysis and optionally uploading results to
   Galaxy. ``illumina_finished_msg.py`` and a Galaxy graphical front end
   both send messages on this queue for starting processing.

-  Called with no queue arguments -- Run individual steps in the
   analysis pipeline. Start multiple servers on distributed machines
   connected with a shared filesystem to allow scaling on a cluster or
   Amazon EC2.

-  Called with ``-q storage`` -- Manage long term storage of larger
   files, like qseq and images.

Specify system specific information in configuration files:

-  ``config/transfer_info.yaml`` -- Configuration on the sequencing
   machine, specifying where to check for new sequencing data.
-  ``config/post_process.yaml`` -- Configuration for analysis and
   storage. This contains links to Galaxy, program locations and
   customization for processing algorithms.
-  ``config/universe_wsgi.ini`` -- Variables used from your Galaxy
   server configuration, including RabbitMQ details for communication
   between the sequencing and analysis machines.

Scripts involved in the processing:

-  ``scripts/bcbio_nextgen.py`` -- Drives the high level analysis of
   sequencing lanes based on information specified through the Galaxy
   LIMS system or in a YAML configuration file. Also produces a PDF
   summary file with statistics on alignments, duplicates, GC
   distribution, quality scores, and other metrics of interest.
-  ``scripts/upload_to_galaxy.py`` -- Handles storing and uploading
   Fastq, alignment, analysis and summary files to Galaxy.

Files generated by the pipeline
~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~

Initial Fastq files (pre-analysis)
**********************************

After basecalling, a number of FastQ files are generated via
solexa\_qseq\_to\_fastq.py script:

::

    1_081227_B45GT6ABXX_1_fastq.txt
    1_081227_B45GT6ABXX_2_fastq.txt
    (...)

The template of them being:

::

    lane_date_fcid_<1|2>_fastq.txt

Where 1\|2 is the forward and reverse read respectively. Only reads that
pass the quality filter (PF) are kept. Quality scores are not converted
to Sanger, but are kept in the produced format (currently Illumina
1.3+). The dots in the fastq files are converted into Ns.

Post-processing generated files
*******************************

Those are distributed in three directories: alignments, barcode and the
top level run directory.

Top level directory
*******************

Those files comprise both plain text files, images and structured data
that can be useful both automatically and in human-readable form to get
different aspects about the sequencing results.

-  ``run_summary.yaml``

Contains a structured view of the run global parameters.

::

    - barcode_id: '2'
      barcode_type: SampleSheet
      lane: '8'
      metrics:
        Aligned: '1204752'
        Pair duplicates: '881332'
        Read length: '101'
        Reads: '11594437'
      request: ''
      researcher: ''
      sample: ''

-  PDF files

Taking a look at an specific sample in a lane, we find the different
sub-components that conform the final summary report
(``*-summary.pdf``). The summary contains GC content, read quality
figures and an insert size histogram.

::

    6_110126_B816J0ABXX_5_1_fastq_txt_quality.pdf
    6_110126_B816J0ABXX_5_2_fastq.txt.segments.hist.pdf
    6_110126_B816J0ABXX_5-sort-dup-gc.pdf
    6_110126_B816J0ABXX_5_1_fastq.txt.segments.hist.pdf
    6_110126_B816J0ABXX_5-sort_1_fastq_qual.pdf
    6_110126_B816J0ABXX_5-sort-dup-insert.pdf
    6_110126_B816J0ABXX_5_2_fastq_txt_quality.pdf
    6_110126_B816J0ABXX_5-sort_2_fastq_qual.pdf
    6_110126_B816J0ABXX_5-sort-summary.pdf

-  BigWig files

Per-sample wigtoBigWig converted files.

-  ``*_metrics files``

They contain plain text values from picard, namely:

::

    1_110126_B816J0ABXX_5-sort-dup.align_metrics
    1_110126_B816J0ABXX_5-sort-dup.gc_metrics
    1_110126_B816J0ABXX_5-sort-dup.dup_metrics
    1_110126_B816J0ABXX_5-sort-dup.insert_metrics

They contain output metrics from Picard, which are parsed later and used
to generate plots. For instance, for ``sort-dup.align_metrics``, the
output from net.sf.picard.analysis.CollectAlignmentSummaryMetrics is
stored.

alignments directory
********************

Contains the results of the alignments for each sample. As we see on the
listing below, lane 1, barcode id 5 has been aligned in `SAM`_ and BAM
formats. For convenience, to facilitate SNP calling, for instance, a
sorted BAM file is also generated.

::

    1_110126_B816J0ABXX_5.sam
    1_110126_B816J0ABXX_5.bam
    1_110126_B816J0ABXX_5-sort.bam
    1_110126_B816J0ABXX_5_1_fastq.bam

_barcode directories
*********************

Those contain fastq files conforming with the naming schema we've seen
before. They are the result of the demultiplexing process, where the
"unmatched" files contain the reads that have not passed the approximate
barcoding matching algorithm:

::

    4_110126_B816J0ABXX_1_1_fastq.txt
    4_110126_B816J0ABXX_5_2_fastq.txt
    4_110126_B816J0ABXX_1_2_fastq.txt
    4_110126_B816J0ABXX_6_1_fastq.txt
    4_110126_B816J0ABXX_5_1_fastq.txt
    4_110126_B816J0ABXX_6_2_fastq.txt

    4_110126_B816J0ABXX_unmatched_1_fastq.txt
    4_110126_B816J0ABXX_unmatched_2_fastq.txt

    4_110126_B816J0ABXX_bc.metrics
    SampleSheet-barcodes.cfg

``*-barcodes.cfg`` contains a simple mapping between barcode id's and
the actual barcode sequence:

::

    3 ATCACGA
    2 ACTTGAA
    9 TAGCTTA
    (...)

The ``_bc.metrics`` file has a plain read distribution for each barcode:

::

    2   11594437
    3   20247932
    9   14390566
    unmatched   908420

Barcodes are added to the 3' end of the first sequence. That way, it
remains platform-independent and can be easily handled downstream. This
`GitHub discussion`_ explains how demultiplexing works. The
demultiplexing is performed by the `barcode\_sort\_trim.py`_ script.

.. _SAM: http://bioinformatics.oxfordjournals.org/content/early/2009/06/08/bioinformatics.btp352.short
.. _GitHub discussion: https://github.com/chapmanb/mgh_projects/commit/3387d82f3496025ad13b69e8d9cbb47cf6ee2af9#nglims_paper/nglims_galaxy.tex-P57
.. _barcode\_sort\_trim.py: https://github.com/brainstorm/bcbb/blob/master/nextgen/scripts/barcode_sort_trim.py
