简介
单细胞领域目前算是进入白热化的阶段,自09年第一篇单细胞文章发表已经度过了10个年头,在通量和成本做到较好的平衡是Department of Genetics, Harvard Medical School,mccarroll实验室15年开发的基于液滴微流控系统生成单细胞Drop-seq技术。目前商业化最成功的10Xgenomics公司开发的单细胞产品也是在此技术的基础上优化而成。HCA人类细胞图谱计划得益与目前这些技术,快速实现了人类多个组织器官的单细胞图谱绘制。本文简要介绍Drop-seq实验的原理及数据分析工具Drop-seq_tools流程要点。
Drop-seq实验
Drop-seq主要核心部件为增压装置及液滴微流控合成装置,原理利用压力将磁珠,细胞,矿物油三种材料,压到微通道里,实现三个物质混成一个大液滴完成单细胞分离标记,完成后续建库。
Drop-seq_tools分析流程
不同的实验方案生成的单细胞数据,分析上整体差别不大,无外乎就是原始数据过滤质控,然后比对到基因组上,根据UMI及细胞barcode实现基因
Drop-seq_tools配合Picardtools,samtools STAR等工具实现对fq bam sam 等文件处理。流程首先将构建参考基因组文件夹,对物种的参考基因组fa及gtf文件做处理。1
2
3
4
5
6
7
8● fasta : The reference sequence of the organism. Needed for most aligners.
● dict : A dictionary file as generated by Picard’s CreateSequenceDictionary . Needed for Picard Tools.
● gtf: The principle file to determine the location of genomic features like genes, transcripts, and exons. Many other metadata files we use derive from this original file. We download our GTF files from ensembl, which has a handy description of the file format here . Ensembl has a huge number of prepared GTF files for a variety of organisms here .
● refFlat: This file contains a subset of the the same information in the GTF file in a different format. Picard tools like the refFlat format, so we require this as well. To make life easy, we provide a program ConvertToRefFlat that can convert files from GTF format to refFlat for you.
● genes.intervals: The genes from the GTF file in interval list format . This file is optional, and useful if you want to go back to your BAM later to see what gene(s) a read aligns to.
● exons.intervals: The exons from the GTF file in interval list format . This file is optional, and useful if you want to go back to your BAM and view what exon(s) a read aligns to.
● rRNA.intervals: The locations of ribosomal RNA in interval list format . This file is optional, but we find it useful to later assess how much of a dropseq library aligns to rRNA.
● reduced.gtf: This file contains a subset of the information in the GTF file, but in a far more human readable format. This file is optional, but can be generated easily by the supplied ReduceGTF program that will take a GTF file as input.
然后用Picard FastqToSam将下机文件转换成sam/bam
随后就可以利用drop-seq 的工具对bam文件做进一步的处理包括:
Example Cell Barcode:1
2
3
4
5
6
7
8
9
10TagBamWithReadSequenceExtended
INPUT=my_unaligned_data.bam
OUTPUT=unaligned_tagged_Cell.bam
SUMMARY=unaligned_tagged_Cellular.bam_summary.txt
BASE_RANGE=112
BASE_QUALITY=10
BARCODED_READ=1
DISCARD_READ=False
TAG_NAME=XC
NUM_BASES_BELOW_QUALITY=1
Example Molecular Barcode:1
2
3
4
5
6
7
8
9
10TagBamWithReadSequenceExtended
INPUT=unaligned_tagged_Cell.bam
OUTPUT=unaligned_tagged_CellMolecular.bam
SUMMARY=unaligned_tagged_Molecular.bam_summary.txt
BASE_RANGE=1320
BASE_QUALITY=10
BARCODED_READ=1
DISCARD_READ=True
TAG_NAME=XM
NUM_BASES_BELOW_QUALITY=1
FilterBAM1
2
3TAG_REJECT=XQ
INPUT=unaligned_tagged_CellMolecular.bam
OUTPUT=unaligned_tagged_filtered.bam
TrimStartingSequence
1 | TrimStartingSequence |
PolyATrimmer1
2
3
4
5
6PolyATrimmer
INPUT=unaligned_tagged_trimmed_smart.bam
OUTPUT=unaligned_mc_tagged_polyA_filtered.bam
OUTPUT_SUMMARY=polyA_trimming_report.txt
MISMATCHES=0
NUM_BASES=6
SamToFastq1
2
3java Xmx4g jar /path/to/picard/picard.jar SamToFastq
INPUT=unaligned_mc_tagged_polyA_filtered.bam
FASTQ=unaligned_mc_tagged_polyA_filtered.fastq
以上完成对比对前的数据做完处理,后续做比对处理。
STAR1
2
3
4/path/to/STAR/STAR
--genomeDir /path/to/STAR_REFERENCE
--readFilesIn unaligned_mc_tagged_polyA_filtered.fastq
--outFileNamePrefix star
SortSam1
2
3
4java Xmx4g jar /path/to/picard/picard.jar SortSam
I=starAligned.out.sam
O=aligned.sorted.bam
SO=queryname
MergeBamAlignment1
2
3
4
5
6
7java Xmx4g jar /path/to/picard/picard.jar MergeBamAlignment
REFERENCE_SEQUENCE=my_fasta.fasta
UNMAPPED_BAM=unaligned_mc_tagged_polyA_filtered.bam
ALIGNED_BAM=aligned.sorted.bam
OUTPUT=merged.bam
INCLUDE_SECONDARY_ALIGNMENTS=false
PAIRED_RUN=false
TagReadWithGeneExon1
2
3
4
5TagReadWithGeneExon
I=merged.bam
O=star_gene_exon_tagged.bam
ANNOTATIONS_FILE=${refFlat}
TAG=GE
DetectBeadSynthesisErrors Detecting and repairing barcode synthesis errors1
2
3
4
5
6
7DetectBeadSynthesisErrors
I=my.bam
O=my_clean.bam
OUTPUT_STATS=my.synthesis_stats.txt
SUMMARY=my.synthesis_stats.summary.txt
NUM_BARCODES= <roughly 2x the number of cells>
PRIMER_SEQUENCE=AAGCAGTGGTATCAACGCAGAGTAC
Digital Gene Expression1
2
3
4
5DigitalExpression
I=out_gene_exon_tagged.bam
O=out_gene_exon_tagged.dge.txt.gz
SUMMARY=out_gene_exon_tagged.dge.summary.txt
NUM_CORE_BARCODES=100
Cell Selection1
2
3
4BAMTagHistogram
I=out_gene_exon_tagged.bam
O=out_cell_readcounts.txt.gz
TAG=XC
后续的数据可以接入Seurat做后续的过滤分群处理。
参考材料:
https://github.com/broadinstitute/Drop-seq
http://mccarrolllab.org/dropseq/
关注微信公众号,第一时间获取更新