实验记录 | 8/5

忙活好现在一堆乱七八糟的事情，到现在已经将近11点了。
那么今天就从现在开始。
首先回顾前两天做的事情。前两天主要就是：
（1）使用CML的其他数据，再次验证在已知标签的情况下，我们能不能将这些样本聚类在一起，也就是说验证sclinager的可行性。==>现在这部分的数据已经在运行，预计8.7结束可以完成这28个样本的运行，并得到阶段性的结果。

• 按照已知类别标签标记的热图，看是否可以聚在一起。

（2）另一方面，想要尝试优化这个流程。最重要的是时间方面的，有无必要砍掉预处理的环节，因为时间消耗最大，并且在有关文献中说其实并无必要。

• 我需要图表方面的论证。将数据的大小以及数据每个环节所用的时间这部分的信息提取出来，绘制成比较直观的图。
• 阅读李恒的那篇文献，理解他的思路，为什么他觉得这一块可以去除？
• 用我们自己的数据进行验证（处理前后的效果有多大变化？）
  （3）了解10X的测序原理，了解我们的数据特征，利用scSplitter这个工具，拆分我们自己的数据（数据很贵，处理的时候要倍加珍惜）。

那么今天我要做的事情是什么？
（1）阅读李恒的那篇文献。我前段时间没有好好看。
（2）根据cellranger的介绍，理解我们现在所用的数据的基本特征。
（3）使用我们昨天拆分的一个细胞的fastq.gz文件，看一下是否可以运行找突变这一流程。

• 我比较担心的是因为我现在还有另外一个流程在使用somatic.pl这个文件，我怕两者之间会相互干扰。
• 我觉得可以先看一下，拆分成了多少个细胞，每一个细胞的数据的大小，估计总的用时。或许是目前更加有效的信息 。

（4）整理出自己的文献阅读的list，目前靶向于variant caller的主要原理，因为我想要优化这个pipeline。

---

ok，现在开始。

1。李恒大佬的文献阅读

这是一篇关于variant calling的review，我对它的兴趣主要有两个方面。
（1）作者是李恒，我挺感兴趣他是如何看待这个问题的。
（2）我想优化QBRC的variant calling 的pipeline，我想看一下，为什么他说可以去除突变前预处理这一环节。
不过这篇文章也有缺点，就是时间有点远（2014）。
话不多说，现在开始。

（11：15）

这张截图，就是我比较关心的问题。我比较好奇的是这个differ by 0.1%的差值是如何计算出来的？

我明显还是读书太少。这张截图就说了不同的caller对于depth的定义是不同的，所以我之前发现了这个问题，产生疑惑的主要原因还是读书太少。

如何理解这里的 variant quality?

2。看一下一共有多少细胞，这些数据有多大。

ls | wc -l #因为目录下全部是文件，没有文件夹所以这条指令能够实现我们的需求
# 25414 
#因为是双端的数据，所以除以2，结果是12707个细胞

统计目录下文件的大小,生成文件large.txt，以后会继续使用。

du -sh *| sort -n >large.txt

输出的文件的内容。

1.0K cell_ACGACTTTTACTTCCT_R1.fastq.gz
1.0K cell_ACGACTTTTACTTCCT_R2.fastq.gz
1.0K cell_ACTAACGAGACTCTGG_R1.fastq.gz
1.0K cell_ACTAACGAGACTCTGG_R2.fastq.gz
1.0K cell_ATTCACTGTGCCAGAC_R1.fastq.gz
1.0K cell_ATTCACTGTGCCAGAC_R2.fastq.gz
1.0K cell_CATGGGTGGTGGTGCA_R2.fastq.gz
1.0K cell_CCGCGGGTGGTGGTGC_R1.fastq.gz
1.0K cell_CCGCGGGTGGTGGTGC_R2.fastq.gz
1.0K cell_CGATTTGTCTGGTTAA_R1.fastq.gz
1.0K cell_CGATTTGTCTGGTTAA_R2.fastq.gz
1.0K cell_CGCATGAGACTCTGGC_R2.fastq.gz
1.0K cell_CTATGGGTGGTGGTGC_R1.fastq.gz
1.0K cell_CTATGGGTGGTGGTGC_R2.fastq.gz
1.0K cell_CTCCGTGGGTGGTGGT_R2.fastq.gz
1.0K cell_CTGGCTTGCTAACTAG_R1.fastq.gz
1.0K cell_CTGGCTTGCTAACTAG_R2.fastq.gz

所以，从某种程度上看这些单个的细胞其实并不是特别的大个。

看一下，解压之后单个细胞的fastq的内容。

(base) [xxzhang@mu02 output]$ head fastq1.fastq

@A00583:531:HW5FTDSXY:2:1668:21658:8093
TTTTTTTTTTGAGACAGTCTCGCTCTTTTACCCAGGCTGGACTGCAGTGGCGCTCTCTCGGCTCACTGCAAGCTCCGCCTCCCAGTCGCCATTCTCCTGCCACAGCCTC
+
FF,FF,FF,FF:FF:,FF:F:F,FFF:F,FFFFF:FFF,FF:::F::::FF,F::::::,,FFFFF:,F,F::,::,FF::FFFF,:FF:F,FFFFF:FFF,FF:FFFF
@A00583:531:HW5FTDSXY:2:1256:16929:28432
TTTTTTTTTTTTTTTTAAACCTTCCTTAAAGATTCTTTGATGCTTTGCTCTATCACTGTAGACCTGGTCTTTTTCCCCCCAGTTTTTTCTTTTTTACATTCTGGGTTG
+
FFFFFFFFFFFFFFFFFFFFF:FF:F:FFFFFF,FFFFFFFF::FFFFFFF,FF:F:,FFF::FFF,FFF:F:FFFFFFFF:F:FF,:::F:,FFFFF:::F,F:F::
@A00583:531:HW5FTDSXY:2:1521:20989:8938
TTTTTTTTTTTTTTTTCGATTCTTTACTTCCTTCATTACACTCTCTATCTATATTAAGTTGTTCTTCTTCTTTTTTAATCTCATCTTTCCACTTTTTATCTATTTTAAT

尝试运行指令，统计运行时间。

perl  /home/xxzhang/workplace/QBRC//somatic.pl  NA  NA  RNA:/home/xxzhang/workplace/QBRC/data/fastq/output/results/cell_ATGAAGCCAAGCCACT_R1.fastq.gz  /home/xxzhang/workplace/QBRC/data/fastq/output/results//cell_ATGAAGCCAAGCCACT_R2.fastq.gz  32  hg38  /home/xxzhang/workplace/QBRC/geneome/hg38/hg38.fasta  /home/xxzhang/miniconda3/bin/java  /home/xxzhang/workplace/QBRC/output_test/  human 1 ./disambiguate_pipeline >pipeline_test.txt

18:33-19:17 1.4M 1个小时。
19:30-

另外我还想看一下这个数据每一个细胞的平均测序深度。

3。下载了Typora，准备了要阅读的文献的list

下午读了lofreq，简直太痛苦了。啥也没读明白。我觉得我挺喜欢读故事的。但是为什么放到文献这边，就有点痛苦的体验了呢？

4。cellranger的网页版的介绍

cellranger这边的信息没啥稀奇的地方，唯一值得注意的是这边的reads的重复率还是挺高的。所以，对于我们而言，可能是在提取细胞数据的过程中，通过调整参数，保留一个UMI的时候删去。也可能是在我们后期mark duplicates的时候删去。

所以，现在今天的任务基本上也是完成了的。虽然我现在腰酸背痛，心情不愉悦。我想到九点的时候就溜。不开心就是不开心，没啥好勉强的。而且现在肚子确实是好饿啊。

现在遗留下来了两个问题。
（1）cutoff值的设定的问题。但是这个cutoff和我们数据的本身的特征相关。是否有必要计算一共平均的coverage的值。我觉得这个平均的coverage应该是针对于我们的exon而言的。
（2）单个细胞运行时间太久，我们是否可以舍去精度去求时间的提高。（两周我们还能够接受，10000个小时实在太可怕了）
（3）有时候会有一种无力感。即使尝试做出一些改变，对于结果的影响也是很小的。
（4）那么这一切的问题，最后的答案回归到了，我们的需求是什么？我们是需要一切尽善尽美，还是虽然这个过程不可控，但是结果能够保证一定的可靠性。
（5）提取细胞信息的这个代码我还没有去看，这也是我心存疑惑的地方？会不会这个环节出了问题呢？
我是不会走的。我要留到九点钟，或许在这一小时内我就有灵感了呢？