混合单细胞rna序列的无基因型解复用
scSplit-jon-xu的Python项目详细描述
SCSPLIT
< H3>混合单细胞RNA SEQ的基因型无解复用，使用隐藏状态模型来识别混合群体内的基因差异样本。它已用于在10x数据集上解复用3到8个样本。
如何安装：
安装python 3.6+
确保可以导入以下python包：
numpy、math、pandas、pickle、pysam、random、scipy、statistics、scikit learn、pyvcf
git克隆https://github.com/jon-xu/scSplit或pip安装scsplit
使用“scsplit<；command>；<；args>；”或“python scsplit<；command>；<；args>；”运行
如何运行工具集：
一。数据质量控制和过滤
a）将目标bam文件（用cb:z:tag标记的条形码）复制到scsplit的同一个文件夹中，只保留带白名单条形码的读取，以减少技术噪音。
b）对来自SCRNA序列的BAM文件进行处理时，应过滤掉以下任何模式的读取：质量低于10，是未映射段，是二级对齐，不是通过过滤器，是PCR或光学复制，或是补充对齐。
示例：samtools view-s-b-q 10-f 3844 original.bam>；target.bam
c）使用samtools中的rmdup、sort、index命令标记BAM文件以进行复制，并对其进行排序和索引
2.呼吁单核苷酸变异
a）在第一步处理后，使用freebayes v1.2从混合样本bam文件中调用snv，设置freebayes的参数，以便不捕获插入和删除（indels）、多核苷酸多态性（mnp）或复杂事件，将最小等位基因计数设置为2，将最小基质量设置为1。
例如：freebayes-f<；reference.fa>；-ixu-c 2-q 1 target.bam>；snv.vcf
这一步可能需要很长时间（如果不使用并行处理，则长达30个小时），gatk或其他snv调用工具也可以工作。用户还可以通过染色体分割bam并分别调用snv。
b）输出VCF文件将进一步过滤，以便仅保留质量分数大于30的SNV。
三。建立等位基因计数矩阵
a）运行“scsplit count”并获取两个.csv文件（“ref_filtered.csv”和“alt_filtered.csv”）作为输出。 输入参数：
-v, --vcf, VCF from mixed BAM
-i, --bam, mixed sample BAM
-b, --barcodes, barcodes whitelist
-r, --ref, Ref count CSV as output
-a, --alt, Alt count CSV as output
e.g.: scSplit count -v mixed_genotype.vcf -i mixed.bam -b barcodes.tsv -r ref_filtered.csv -a alt_filtered.csv
b）此步骤消耗内存，所需的ram高度依赖于上一步的snv数量和单元数量。作为指导，一个包含60000个snv和10000个单元格的矩阵可能需要超过30gb的ram才能运行，请允许足够的ram资源来运行脚本。
c）可使用两种滤波方法对矩阵进行滤波，以提高预测精度：
重复序列区域（repeatMasker）
HG19:http://www.repeatmasker.org/genomes/hg19/RepeatMasker-rm405-db20140131/hg19.fa.out.gz
HG38:http://www.repeatmasker.org/genomes/hg38/RepeatMasker-rm405-db20140131/hg38.fa.out.gz
若要处理下载的文件，请获取文件的染色体、开始列和结束列，并使用这些黑名单区域筛选被调用的snv，并获取不属于重复区域的snv列表（格式为：“1:10177”）。hg19和hg38的已处理重复区域可在以下位置找到：http://data.genomicsresearch.org/Projects/scSplit/RepeatRegions
普通单核苷酸多态性（例如，1000基因组计划中的人类普通单核苷酸多态性）
HG19:ftp://ftp.1000genomes.ebi.ac.uk/vol1/ftp/release/20130502/
要处理常见snp的基因型文件，可以下载每个染色体文件，并使用bcftools或整个基因组文件将它们连接起来，取vcf文件的前两列并用冒号替换制表符，使每一行都是一个snv，例如“1:10177”。hg19和hg38的已处理通用snv可在以下位置找到：http://data.genomicsresearch.org/Projects/scSplit/CommonSNVs
然后使用公共snv列表筛选在最后一步（“ref_filtered.csv”和“alt_filtered.csv”）中生成的矩阵，并将它们用作参考矩阵和备用矩阵作为scsplit运行的输入。
四。解复用并生成ALT P/A矩阵< EH3>
a）使用两个生成的等位基因计数矩阵文件将细胞分解成不同的样本。doublet sample不会每次都有相同的sample id，这将在日志文件中明确指出
b）这一步也消耗内存，所需的ram高度依赖于上一步的snv数量和单元数量。作为指导，一个包含60000个snv和10000个单元格的矩阵可能需要50gb以上的ram才能运行，请允许足够的ram资源来运行脚本。
c）运行“scsplit run”，输入参数：
-r, --ref, Ref count CSV as input
-a, --alt, Alt count CSV as input
-n, --num, Number of mixed samples
-v, --vcf, individual genotypes to check distinguishing variants against (optional)
e.g.: scSplit run -r ref_filtered.csv -a alt_filtered.csv -n 8
d）将生成以下文件：
"scSplit_doublet.txt": indicating which cluster is doublet state
"scSplit_barcodes_{n}.csv": N+1 indicating barcodes assigned to each of the N+1 samples (including doublet state)
"scSplit_dist_alleles.txt": the distinguishing alleles that can be used to genotype and assign sample to clusters
"scSplit_dist_matrix.csv": the ALT allele Presence/Absence (P/A) matrix on distinguishing variants for all samples as a reference in assigning sample to clusters
"scSplit_PA_matrix.csv": the full ALT allele Presence/Absence (P/A) matrix for all samples
"scSplit.model", a python pickle dump containing the final allele fraction model (model.model_MAF), and the probability of each cell belonging to each sample (model.P_s_c)
"scSplit.log" log file containing information for current run, iterations, and final Maximum Likelihood and doublet sample
5个。根据分割结果生成样本基因型
a）运行“scsplit基因型”，输入参数：
-r, --ref, Ref count CSV as output
-a, --alt, Alt count CSV as output
-p, --psc, generated P(S|C)
e.g.: scSplit genotype -r ref_filtered.csv -a alt_filtered.csv -p scSplit_P_s_c.csv
b）将为所有样本模型的对数转换基因型likelihood生成VCF文件（“scsplit.vcf”）。
欢迎加入QQ群-->： 979659372
