Bulk / Smart-seq2数据分析一般流程

Bulk RNA-seq（RNA-Seq of bulk samples）是一种RNA-Seq技术应用，它是通过将整个组织或细胞群体的RNA提取并混合，进行高通量测序来分析基因表达的技术。

Smart-seq2是一种单细胞转录组测序技术，用于对单个细胞中的RNA进行测序，并可以对单个细胞的基因表达进行分析。在流程上，Smart-seq与Bulk RNA-seq的一般流程基本相同。

以下是Bulk/Smart-seq2数据分析的主要流程：

1. 过滤：首先需要对原始的测序数据进行质控和过滤，包括去除低质量的reads、去除含有接头污染的reads、去除低复杂度的reads等。

# 过滤
fastp -i input.fastq -o output.fastq -h report.html -j report.json

2. 质控：使用fastqc对过滤后的数据进行质控，检查过滤效果。

#质控
fastqc -t 20 -o fastqc ./all_data/*fq.gz

3. 读取比对：将过滤后的reads比对到参考基因组上，可以使用常见的比对工具，如STAR、Hisat2等。比对后需要对比对结果进行质量评估，例如计算比对率、可重复性、覆盖率等指标。

# 比对到参考基因组
 hisat2 -p 20 \
-x /path to ref/ref_Mmul_10_gtf \
-1 ./fastp/Data1.fq.gz \
-2 ./fastp/Data2.fq.gz \
-S ./bam/data.sam > ./bam/data_hisat2.log 2>&1

4. 转录本重建和定量：使用转录本定量工具，如StringTie、Cufflinks等，对比对后的reads进行转录本定量。转录本定量结果可以用于后续的差异表达分析和聚类分析。

# 转录本重建和定量
stringtie -p 8 -G /path/to/annotation.gtf -o output.gtf -l output input.bam

5. 差异表达分析：对单细胞测序数据进行差异表达分析，可以使用常见的R软件包，如DESeq2、edgeR等。通过比较不同细胞之间的基因表达差异，可以发现关键差异表达基因和细胞亚型等信息。

# 差异表达分析
library(DESeq2)
counts <- read.table("counts.txt", header = T, row.names = 1)
metadata <- read.table("metadata.txt", header = T, row.names = 1)
dds <- DESeqDataSetFromMatrix(countData = counts, colData = metadata, design = ~ condition)
dds <- DESeq(dds)
res <- results(dds)

6. 功能注释和富集分析：对差异表达基因进行功能注释和富集分析，以帮助理解差异表达基因的功能和调控网络，从而获得对研究样品的更深入的理解。（以KEGG为例）

# 使用biomaRt包从Ensembl数据库中获取基因注释信息
library(biomaRt)
ensembl <- useMart("ensembl")
mart <- useDataset("mmusculus_gene_ensembl", mart = ensembl)
genes <- c("ENSMUSG00000026561", "ENSMUSG00000044791", "ENSMUSG00000021584")
annotations <- getBM(attributes = c("ensembl_gene_id", "entrezgene_id", "description"),
                      filters = "ensembl_gene_id",
                      values = genes,
                      mart = mart)
# 使用KEGGREST包从KEGG数据库中获取通路信息
library(KEGGREST)
pathways <- keggList("pathway", "mmu")
pathway_genes <- list()
for (pathway in pathways) {
  genes_in_pathway <- keggGet(pathway, "genes")
  pathway_genes[[pathway]] <- genes_in_pathway
}

7. 数据可视化：最后，将聚类分析的结果进行可视化，可以使用常见的数据可视化工具，如ggplot2、Plotly等，来帮助用户更好地理解数据。

# 数据可视化
library(ggplot2)
ggplot(data = res, aes(x = log2FoldChange, y = -log10(pvalue))) +
  geom_point() +
  theme_classic() +
  xlab("log2FoldChange") +
  ylab("-log10(pvalue)")