windows ubuntu 子系统：肿瘤全外篇,bam质控

各个环节的质控，
raw和clean都要质控，
比对的各环节的bam文件都要质控，
使用qulima对wes的比对bam文件总结测序深度及覆盖率。

samtools flagstat L1_recalibrated_reads.bam
该命令将输出 BAM 文件的一些统计信息，包括总读取数、比对上参考序列的读取数、比对到不同位置的读取数等。

#结果可如下。

L1_recalibrated_reads.bam 的统计信息如下：
总读取数：103,094,432
比对上参考序列的读取数：103,028,917 (占总读取数的 99.94%)
次要比对的 reads 数：0
补充比对的 reads 数：674,520
重复 reads 数：22,411,852
成对测序的 reads 数：102,419,912
测序的 read1 数：51,209,956
测序的 read2 数：51,209,956
正确成对匹配的 reads 数：101,697,064 (占成对测序的 reads 的 99.29%)
自身及其 mate 均比对到参考序列的 reads 数：102,306,392
单独出现的 reads 数：48,005 (占总读取数的 0.05%)
与不同染色体的 mate 均比对的 reads 数：400,816
映射到不同染色体且 mapQ 值大于等于 5 的 reads 数：304,376

运行以下命令可以计算 L1.bam 中的总行数（即记录数），从而得知该 BAM 文件中包含多少条比对信息：samtools view 949743-T_L2_1.bam | wc -l

#获取全外bed文件

CCDS官网
进入官网后进入其ftp服务器

cat CCDS.20221027.txt | perl -alne '{/

$$(.*?)$$ /;next unless $1;$gene=$F[2];$exons=$1;$exons=~s/\s//g;$exons=~s/-/\t/g;print "$F[0]\t$_\t$gene" foreach split/,/,$exons;}'|sort -u |bedtools sort -i |awk '{print "chr"$0"\t0\t+"}'  > hg38.exon.bed
这条命令的作用是将 CCDS（Consensus CDS）数据中的 exons 信息提取出来，生成一个 BED 文件 hg38.exon.bed。具体实现步骤如下：使用 cat 命令将 CCDS.20221027.txt 文件的内容输出到标准输出。
使用 perl 命令解析每一行，并通过正则表达式提取出 exons 信息。如果该行不包含 exons 信息，则跳过。
将提取到的 exons 信息进行格式化，并使用 split 函数将其拆分成多个 exon。对于每个 exon，输出其所在的染色体、起始位置、终止位置和所属基因。
使用 sort 命令将输出结果按照染色体、起始位置和终止位置排序。
使用 awk 命令将排序后的结果转换为 BED 格式，并指定其 score 和 strand 信息，最终将结果输出到 hg38.exon.bed 文件中。
这个 hg38.exon.bed 文件可以用于基因组注释和区域相关的分析。

samtools view L1_recalibrated_reads.bam | less -S
这条命令使用 samtools view 命令来查看 949743-T_L2_1_recalibrated_reads.bam 这个 BAM 文件的内容，并通过管道将输出传递给 less -S 命令进行分页查看。
samtools view 命令用于从 BAM 文件中读取比对信息，并以文本格式输出。| 符号表示将前一个命令的输出作为后一个命令的输入进行处理。
less 命令是一个分页查看器，可以按需滚动查看文件的内容。-S 参数用于禁用行内过长时的折行显示，保持每行内容在屏幕上的可见性。
因此，执行该命令后，将能够使用 less 分页查看 L1_recalibrated_reads.bam 文件中的比对信息。您可以使用方向键（上下左右）和 Page Up/Page Down 键来浏览文件内容，并使用 q 键退出 less 查看器。

# 1. 创建输出目录
mkdir -p qc_results

#安装qualimap

qualimap bamqc \
    -bam L1.bam \
    -outdir qc_results \
    -c \
    --java-mem-size=4G \
    --feature-file /mnt/h/db/hg38.bed/hg38.exon.bed \
    -nt 4

qualimap bamqc: 这是运行 Qualimap 工具中的 bamqc 模块的命令，用于评估 BAM 文件的质量。

-bam L1.bam: -bam 参数指定输入的 BAM 文件，这里使用的是 949743-T_L2_1.bam 文件。

-outdir qc_results: -outdir 参数指定输出结果的目录，这里结果将保存在名为 qc_results 的目录中。

-c: -c 参数表示生成覆盖度报告。

--java-mem-size=4G: --java-mem-size 参数指定分配给 Java 虚拟机的内存大小，这里设置为 4GB。

--feature-file /mnt/h/db/hg38.bed/hg38.exon.bed: --feature-file 参数指定感兴趣的区域文件，这里使用的是一个 BED 格式的文件，其中包含了人类基因组 hg38 版本的外显子区域信息。

-nt 4: -nt 参数指定并行运行的线程数，这里设置为 4 个线程。

出来以下结果，有些难懂。

可用multiqc整理一下就好看多了。