Introduction

rl/hifiampliconanalysis is a bioinformatics pipeline that takes amplicon data generated from PacBio HiFi reads, runs quality control on the reads, and aligns the reads to a reference gene and/or cDNA sequence. Reads from different genes should be separated, as the pipeline processes amplicons from a single gene, but can parallelize across multiple different genes. These alignments are then used to cluster reads based on sequence similarity into distinct sub-groups representing the different amplicon species present. Variants are called and a consensus sequence is output (along with a VCF of high confidence variants) for each cluster that pass stringent filtering criteria. A linear plot is output, showing each different species of amplicon present along with variants and gene structure.

1. Read QC (Nanoplot)
2. Convert fastq to fasta (seqtk)
3. Screen and remove contaminants (FCS-adaptor)
4. Alignment (Minimap2, samtools)
5. Filter reads without primers (Biopython)
6. Cluster sequences (Decipher)
7. Align clusters (Minimap2, samtools)
8. Tag bams (pysame)
9. Split bams by cluster (samtools)
10. Call variants (bcftools)
11. Call consensus sequences (bcftools)
12. Find the longest ORF (orfipy)
13. Plot amplicon species (ggtranscript)

Usage

[!NOTE] If you are new to Nextflow and nf-core, please refer to this page on how to set-up Nextflow. Make sure to test your setup with -profile test before running the workflow on actual data.

First, prepare a samplesheet in the assets folder with your input data that looks as follows, ensuring that the necessary files are placed in a data/ folder:

samplesheet.csv:

csv sample,fastq,fasta,primer1,primer2,bed sample1,data/sample1.fastq.gz,data/reference.fasta,ACTG...,ACTG...,data/reference.bed

Each row represents a fastq file (PacBio HiFi), a gene reference file, the sequences of the gene-specific primers used for amplification, and a bed file that has the chromosome and positions of the gene features of interest (introns/exons/UTRs, etc...). The bed file should be derived from the same gene/cDNA sequence as the reference fasta and can be obtained by navigating to the UCSC genome browser, searching for the gene of interest, clicking the Tools dropdown, going to the Table Browser, then selecting the region of interest along with the knownGene table, selecting BED as the output format, select Get output, then selecting and downloading each feature you are interested in, followed by collating all features into a single BED file.

Now, you can run the pipeline using:

bash nextflow run rlinder02/rl-hifiampliconanalysis \ -r v1.0.0 -profile <docker/singularity/.../institute> \ --input assets/samplesheet.csv \ --outdir <OUTDIR>

[!WARNING] Please provide pipeline parameters via the CLI or Nextflow -params-file option. Custom config files including those provided by the -c Nextflow option can be used to provide any configuration except for parameters; see docs.

Credits

rlinder02/rl-hifiampliconanalysis was originally written by Rob Linder.

Contributions and Support

If you would like to contribute to this pipeline, please see the contributing guidelines.

Citations

An extensive list of references for the tools used by the pipeline can be found in the CITATIONS.md file.

Bo Hu, Jinpu Jin, An-Yuan Guo, He Zhang, Jingchu Luo and Ge Gao. GSDS 2.0: an upgraded gene feature visualization server. Bioinformatics 2015, 31(8):1296-1297. DOI: https://doi.org/10.1093/bioinformatics/btu817.

This pipeline uses code and infrastructure developed and maintained by the nf-core community, reused here under the MIT license.

The nf-core framework for community-curated bioinformatics pipelines.

Philip Ewels, Alexander Peltzer, Sven Fillinger, Harshil Patel, Johannes Alneberg, Andreas Wilm, Maxime Ulysse Garcia, Paolo Di Tommaso & Sven Nahnsen.

Nat Biotechnol. 2020 Feb 13. doi: 10.1038/s41587-020-0439-x.