这次用到的数据集是GSE274995，里面包含了3个样本的头颈部鳞癌细胞系(Cal27细胞)数据。数据上传者采用的上游分析软件：
使用了Trimmomatic进行了数据质控；Bowtie2进行了比对分析；SAMtools处理比对后的BAM文件；MACS2用于Peaks Calling；bamCoverage转换为二进制的BigWig文件; 使用了hg19作为参考基因组，现在基本上首先使用GRCh38。再来看一下曾老师推文和B站视频中的分析流程，其中数据数据过滤环节使用了fastp或者trim-galore软件(替代Trimmomatic)，其他的基本差不多。
总体的步骤可以参考既往mRNA/miRNA的步骤：环境部署——数据下载——查看数据(非质控)——数据质控清洗——数据比对——数据定量，基本逃不出这个框架，只是针对不同类型的数据稍有不同。
本次分析步骤包括：环境部署——数据下载——查看数据(非过滤)
分析步骤
1.环境部署/软件安装：
尝试使用ARM架构(M1/M2芯片)去安装fastqc trim-galore hisat2 subread multiqc samtools salmon fastp，但这些软件中有几个是不兼容的。所以需要改回原来的x86_64架构(Intel芯片)，如果非mac/M1/M2的不需要用这种方式。
# macbook
CONDA_SUBDIR=osx-64 conda create -n chipseq_x86_64 python=3.9
conda activate chipseq_x86_64
# 如果是在服务器中
conda create -n chipseq python=3.9
conda activate chipseq
如果安装软件很慢的话可以先安装mamba再安装软件，比如：只需要把conda修改成mamba即可
conda install -c conda-forge mamba
mamba install -y sra-tools trim-galore samtools deeptools homer meme macs2 bowtie bowtie2
2.安装软件
conda install -y sra-tools trim-galore samtools deeptools
conda install -y homer meme macs2 bowtie bowtie2
3.下载数据
# 建议在自己的服务器或者本地先建好相应的文件夹
# 然后需要先cd到目标文件夹下
# 不建议直接复制，因为可能会有不知名的空格hhh
nohup cat SRR_Acc_List.txt | while read id; do prefetch $id; done > download.log 2>&1 &
nohup的作用是使后续命令在后台运行，即使关闭中断或断开会话任务也会继续进行。 cat可以读取文件的内容并将其逐行输出 “;”和“|”区别是前者可以分隔独立命令，每条命令按顺序执行，前后无依赖关系。后者是链接命令，前一个命令的输出是后一个命令的输入，命令之间存在数据流依赖。 这里为什么可以用id去作为读取对象，是因为read的命令设计会从输入流中读取每一行的数据，并存储到指定的变量中，id是自定义的变量名，假设换成line或者其他都可以。 ">"的作用是输入重定向，也就是说如果在download.log前面使用了比如“|”，这个log文件会被当做命令。 “2>&1”的含义是2(表示标准错误，stderr，文件描述2)，>(表示重定向)，&1(表示重定向的目标是标准输出，stdout, 文件描述1)。 最后的“&”是将整个命令放在后台运行。
4.sra文件转化为fastq
先看看样本的单双链信息，点击具体样本从这里可以看到样本采用的测序方式是什么, Layout: SINGLE,说明是单端测序数据
此外在这里稍提一个概念，有时候我们在让公司帮忙做测序的时候，公司人员交流会提到“测多少个G”这个概念。这里的测多少个G是表示Gigabases(表示Giga碱基，10亿碱基)。
我们可以看到下面这个页面上最上面的这句话：1 ILLUMINA (Illumina NovaSeq 6000) run: 42.3M spots, 3.2G bases, 1.2Gb downloads。
ILLUMINA (Illumina NovaSeq 6000): ILLUMINA：指使用的测序技术平台，Illumina 是一种高通量测序技术。Illumina NovaSeq 6000：是 Illumina 平台的一种具体型号。 run:表示这是一次测序实验（run）的统计数据。 42.3M spots:42.3M：表示测序的 reads（读取数），单位是百万（M，即 10^6）。spots：指的是测序数据中的原始数据点，通常等同于测序读取数（reads），或者是指成对读取的原始数据。 3.2G bases:3.2G：表示测序生成的碱基数，总量是 3.2 亿个碱基（G = 10^9）。bases：指测序生成的碱基序列。这里G bases的简写也是 Gb，要与下面的Gigabytes作区分，这里表示碱基数。 1.2Gb downloads:1.2Gb：表示这次测序运行生成的文件大小约为 1.2 吉字节（Gb = gigabytes）。downloads：指该测序运行的数据文件总大小为 1.2GB，可以通过 SRA 下载，这里的Gb表示存储量大小。 一般来说，测序数目(42.310^6) * 测序长度(?这里是未知) = 测序量(3.210^9 bases)，我们反推一下测序长度，当然这个数据是单端测序数据，如果是双端我们还需要除以2。 我们用一种不太智能的方法验证一下是不是75bp，我们进入如下界面，点开reads。
数一数确实是75个碱基。
接下来正式开始进行fastq转换
# 由于笔者这里存储的有点乱，所以会比较复杂
# 使用者在自行使用的时候，如果觉得设置路径很麻烦就把文件全放一个文件夹下面，先学会再说
# 定义变量，指定工作目录路径
# path路径需要通过pwd确定
# 可以选择写一个shell脚本
# 把下面代码输入进去
#!/bin/bash
# 也可以直接复制这些代码进行运行
path="/home/lm/Z_Projects/chipseq"
list_file="${path}/SRR_Acc_List.txt" # 列表文件路径
sra_dir="${path}/sra/sra" # 输入的 .sra 文件所在路径
fastq_dir="${path}/fastq" # 输出的 .fastq 文件存储路径
# 创建 fastq 输出目录
mkdir -p "${fastq_dir}"
# 读取列表文件并处理
cat "${list_file}" | while read id; do
# 检查 .sra 文件是否存在
sra_file="${sra_dir}/${id}.sra"
if [ -f "${sra_file}" ]; then
echo "Processing ${id}..."
# 使用 fastq-dump 处理 .sra 文件，输出到 fastq 文件夹
fastq-dump --gzip -O "${fastq_dir}" "${sra_file}"
echo "Completed ${id}."
else
echo "File ${sra_file} not found, skipping..."
fi
done
读取列表文件并循环处理：cat "{list_file} 变量指定的列表文件内容，通常是一个包含若干 SRR ID（例如 SRR123456）的文本文件，每一行是一个 ID。|: 管道符，将 cat 的输出作为后续 while 循环的输入。while read id: 按行读取列表文件中的内容，并将每一行存储到变量 id 中，逐一处理。do: 标志循环的开始，后续是循环体的内容。 定义 .sra 文件的路径：{id}: 当前循环中读取的 ID 值，例如 SRR123456。"{id}.sra": 拼接出完整的 .sra 文件路径，例如 /home/user/chipseq/sra/sra/SRR123456.sra。 检查 .sra 文件是否存在：[ -f "{sra_file}": 表示要检查的文件路径。 打印处理消息：echo: 打印消息到终端，告知当前正在处理的 ID。 使用 fastq-dump 工具处理文件：fastq-dump: 是一个常用的 SRA 数据转换工具，用于将 .sra 格式转换为 .fastq 或 .fastq.gz 格式。--gzip: 转换后的 .fastq 文件会直接压缩为 .fastq.gz 格式，节省存储空间。-O "{fastq_dir} 目录中。"${sra_file}": 要处理的 .sra 文件路径。 打印完成消息：echo: 打印消息到终端，告知用户当前 ID 的处理已经完成。文件不存在时的提示：else: 如果条件 if 不成立（即文件不存在），执行该分支的代码。echo: 打印警告消息到终端，告知用户某个 ID 的 .sra 文件未找到，并跳过该文件的处理。 循环结束标志：fi: 结束 if 条件判断块。done: 表示 while 循环的结束。
4.查看数据(非过滤)
# cd到fastq文件夹中
# 不整合，输出每一个样本的质控报告
# fastqc -t [线程数] [存储路径] [文件来源路径]
fastqc -t 6 -o ./ ./SRR*.fastq.gz
不整合，输出每一个样本的质控报告
fastqc运行 FastQC 软件; -t 6 运行6个线程，加快速度；./ 存在当前文件夹中； ./SRR*.fastq.gz 读取fastq.gz文件可以打开html文件看一看，具体的质控报告细节可以看一下既往的转录组上游分析流程~ 也可以看看其他UP主的推文。
# 整合，输出一个总的质控报告
# 请注意，这里是先要进行单一质控，然后根据单一质控得到的zip文件结果进行整合！
multiqc ./*.zip -o ./
整合，输出一个总的质控报告。这里是先要进行单一质控，然后根据单一质控得到的zip文件结果进行整合！
既往推文
转录组上游分析流程(一）：转录组上游分析流程(一） 转录组上游分析流程(二）：转录组上游分析流程(二） 转录组上游分析流程(三）：转录组上游分析流程(三） 转录组上游分析流程(四）：转录组上游分析流程(四）
参考资料
生信技能树：https://mp.weixin.qq.com/s/KOFd60USQdKY2C4Sc4OIvA 生信学习日记：https://mp.weixin.qq.com/s/z-qbVZqo-cKT3Tfv4LT6OQ Tobin笔记：https://mp.weixin.qq.com/s/Uf8GgK4kZAWbRKKoh6IiyA
致谢：感谢曾老师以及生信技能树团队全体成员。
注：若对内容有疑惑或者有发现明确错误的朋友，请联系后台(欢迎交流)。更多内容可关注公众号：生信方舟
- END -