玉米中的元基因调控网络突出了功能上干系的调控相互作用\mo.20a3.R ...

海哥  金牌会员 | 2024-12-30 16:09:41 | 显示全部楼层 | 阅读模式
打印 上一主题 下一主题

主题 888|帖子 888|积分 2664

这个脚本的重要作用是处置惩罚和分析RNA-seq数据,进行样本的聚类分析、降维可视化,并且将结果生存下来。详细来说,它包括以下几个关键步骤:
1. 读取RNA-seq数据



  • 利用rnaseq_cpm_raw()和rnaseq_cpm()两个函数分别加载原始的RNA-seq数据和颠末处置惩罚后的RNA-seq数据。
  • 数据包含基因表达矩阵(tm)、样本元数据(th)以及一些基因注释信息(tl)。
2. 样本信息的加载与归并



  • 从samples.xlsx文件中读取样本信息,并与RNA-seq数据中的样本元数据(th)归并。
  • 将样本的基因型、处置惩罚范例、冷耐心等信息拼接成一个新的标签(lab),用于后续分析和可视化。
  • 对RNA-seq数据中的表达量进行变更(利用反双曲正弦变更 asinh(CPM)),以平滑表达数据并减少大数值的影响。
3. 聚类分析(层次聚类)



  • 利用层次聚类(hclust)分析样本之间的相似性,并根据干系性(皮尔逊或斯皮尔曼干系系数)计算基因表达的相似度。
  • 生成聚类结果的可视化图(PDF格式),通过不同的干系性度量(如皮尔逊或斯皮尔曼)来比较样本的聚类情况。
4. 降维分析与可视化



  • 通过**PCA(主成分分析)t-SNE(t分布随机邻域嵌入)**方法将高维的基因表达数据降维到2D空间,帮助可视化数据中的样本之间的差别或聚类。
  • 利用plot_pca()和plot_tsne()生成PCA和t-SNE的可视化图,并生存为PDF文件。
  • 在降维图中,样本的颜色和形状可以代表不同的处置惩罚组或冷耐心状态,帮助辨认不同的样本群体。
5. 样本数据的清算与生存



  • 利用complete_sample_list()函数对样本数据进行清算,去除不完整或无效的样本。
  • 将终极的清算过的样本元数据(th2)生存为TSV格式的文件(01.meta.tsv),以便后续利用或备份。
6. 文件生存



  • 通过ggsave()将各类分析结果(如聚类图、PCA图、t-SNE图等)生存为PDF文件,这些文件可以用于报告或进一步分析。

总结:

该脚本的焦点目标是处置惩罚RNA-seq数据并进行多种数据分析和可视化,包括:


  • 数据预处置惩罚:归并样本信息、处置惩罚缺失值和进行数据变更。
  • 聚类分析:通过层次聚类评估样本间的相似性。
  • 降维与可视化:通过PCA和t-SNE展示样本之间的关系。
  • 结果生存:将聚类和降维的可视化结果生存为文件,并将清算后的样本元数据生存为TSV文件。
这些步骤为后续的生物学分析(如基因表达模式的分析、差别分析等)奠定了基础。
第一部门:设置工作目次并加载数据

  1. genome = 'Osativa'
  2. t_cfg = read_projects(genome)
  3. t_cfg %>% print(n=50)
  4. yid = 'mo20a3'
  5. dirw = glue("{dird}/11_qc/{genome}/{yid}")
  6. if(!dir.exists(dirw)) dir.create(dirw)
复制代码
解析

  • genome = 'Osativa'
    设置了基因组的名称,这里选择的是‘Osativa’(应该是指水稻基因组)。
  • t_cfg = read_projects(genome)
    read_projects()函数用于读取与指定基因组干系的项目配置文件,并将其存储在t_cfg中。通常这里会包含数据的元信息(例如样本、条件等)。
  • t_cfg %>% print(n=50)
    打印t_cfg的前50行信息,通常用于检察该基因组下全部项目标概要。
  • yid = 'mo20a3'
    设置实验的样本ID(yid),代表一个特定的实验或处置惩罚。
  • dirw = glue("{dird}/11_qc/{genome}/{yid}")
    利用glue()函数生成一个文件路径,该路径将存放QC(质量控制)结果。dird是之前定义的一个变量,大概表示根目次或数据文件夹。
  • if(!dir.exists(dirw)) dir.create(dirw)
    查抄是否存在上面创建的路径dirw,如果没有,则创建该目次。
这一部门的重要功能是初始化工作情况并加载数据配置文件,之后为后续的分析步骤设置目次结构。
第二部门:读取原始数据并归并样本信息

  1. #{{{ raw: read in
  2. res = rnaseq_cpm_raw(yid, genome)
  3. th = res$th; tm = res$tm; tl = res$tl; th_m = res$th_m; tm_m = res$tm_m
  4. th2 = read_xlsx(glue("~/projects/cold/data/samples.xlsx")) %>%
  5.     replace_na(list(Name='?', ColdTolerant='?'))
  6. th = res$th %>% left_join(th2, by="Genotype") %>%
  7.     mutate(lab = str_c(Genotype, Name, Treatment, ColdTolerant, sep='_'))
  8. tm = res$tm %>% filter(SampleID %in% th$SampleID) %>%
  9.     mutate(value=asinh(CPM))
  10. #}}}
复制代码
解析

  • res = rnaseq_cpm_raw(yid, genome)
    通过调用rnaseq_cpm_raw()函数,加载原始的RNA-seq数据。返回的res包含了几个重要数据框:

    • th:包含样本的元数据(如样本ID、处置惩罚组等)。
    • tm:包含每个样本的基因表达数据(每个基因的CPM值)。
    • tl:包含基因注释信息。
    • th_m、tm_m:是针对特定的处置惩罚或条件进行过计算的干系数据。

  • th2 = read_xlsx(glue("~/projects/cold/data/samples.xlsx"))
    从Excel文件~/projects/cold/data/samples.xlsx读取样本信息表格,th2是读取的表格数据。
  • replace_na(list(Name='?', ColdTolerant='?'))
    利用replace_na()函数,将缺失的样本名称(Name)和冷耐心标识(ColdTolerant)添补为问号?。
  • th = res$th %>% left_join(th2, by="Genotype")
    利用left_join()将res$th中的元数据和从Excel读取的th2表格数据进行归并,归并的依据是Genotype字段(基因型)。如允许以把样本的基因型信息与额外的形貌性数据联合起来。
  • mutate(lab = str_c(Genotype, Name, Treatment, ColdTolerant, sep='_'))
    创建一个新的列lab,将Genotype、Name、Treatment和ColdTolerant列的内容拼接在一起,并用下划线分隔。这列通常作为标签,用于标识样本。
  • tm = res$tm %>% filter(SampleID %in% th$SampleID)
    在表达数据res$tm中筛选出与归并后的样本信息th中SampleID匹配的行,确保只保留在th中列出的样本。
  • mutate(value=asinh(CPM))
    对表达矩阵中的CPM值应用反双曲正弦(asinh)变更。这是常见的对基因表达数据进行平滑处置惩罚的方式,尤其是在CPM值有较大差别时。
这一部门的关键任务是:


  • 加载和归并样本元数据与表达数据。
  • 对表达数据进行须要的预处置惩罚(如取对数变更)。
第三部门:聚类分析与可视化(hclust & tSNE)

  1. #{{{ raw: hclust & tSNE
  2. p1 = plot_hclust(tm, th, pct.exp = .7, cor.opt = 'pearson', var.col = 'Treatment', expand.x = .4)
  3. ggsave(file.path(dirw, '11.hclust.p.pdf'), p1, width = 6, height = 8)
  4. p1 = plot_hclust(tm, th, pct.exp = .7, cor.opt = 'spearman', var.col = 'Treatment', expand.x = .4)
  5. ggsave(file.path(dirw, '11.hclust.s.pdf'), p1, width = 6, height = 8)
  6. p2 = plot_pca(tm, th, pct.exp = .7, pca.center = T, pca.scale = F,
  7.     var.shape = 'ColdTolerant', var.col = 'Treatment', var.lab = 'Name',
  8.     legend.pos = 'bottom.left', legend.dir = 'v', pal.col = 'aaas')
  9. ggsave(file.path(dirw, '11.pca.pdf'), p2, width = 6, height = 6)
  10. p3 = plot_tsne(tm, th, pct.exp = .8, perp = 3, iter = 1200, seed = 2,
  11.     var.shape = 'ColdTolerant', var.col = 'Treatment', var.lab = 'Name',
  12.     legend.pos = 'top.right', legend.dir = 'v', pal.col = 'aaas')
  13. ggsave(file.path(dirw, '11.tsne.pdf'), p3, width = 6, height = 6)
  14. #}}}
复制代码
解析

  • 层次聚类(hclust)
    1. p1 = plot_hclust(tm, th, pct.exp = .7, cor.opt = 'pearson', var.col = 'Treatment', expand.x = .4)
    2. ggsave(file.path(dirw, '11.hclust.p.pdf'), p1, width = 6, height = 8)
    复制代码

    • plot_hclust()函数用于进行层次聚类分析。
    • tm是基因表达矩阵,th是样本信息。
    • pct.exp = .7 表示只利用表达量超过70%的基因进行聚类分析。
    • cor.opt = 'pearson' 选择利用皮尔逊干系系数来计算基因间的相似性。
    • var.col = 'Treatment' 用样本的“Treatment”字段来为不同样本上色,以便在聚类结果中区分。
    • expand.x = .4 控制图像的x轴扩展程度。
    • 最后利用ggsave()将结果生存为PDF文件。
    之后的部门:
    1. p1 = plot_hclust(tm, th, pct.exp = .7, cor.opt = 'spearman', var.col = 'Treatment', expand.x = .4)
    2. ggsave(file.path(dirw, '11.hclust.s.pdf'), p1, width = 6, height = 8)
    复制代码
    这段代码执行了雷同的聚类分析,不过利用了斯皮尔曼干系系数(cor.opt = 'spearman')来代替皮尔逊干系系数。

  • 主成分分析(PCA)
    1. p2 = plot_pca(tm, th, pct.exp = .7, pca.center = T, pca.scale = F,
    2.     var.shape = 'ColdTolerant', var.col = 'Treatment', var.lab = 'Name',
    3.     legend.pos = 'bottom.left', legend.dir = 'v', pal.col = 'aaas')
    4. ggsave(file.path(dirw, '11.pca.pdf'), p2, width = 6, height = 6)
    复制代码

    • plot_pca()函数执行PCA分析,减少数据维度并可视化。
    • pct.exp = .7 表示利用70%的基因表达数据进行分析。
    • pca.center = T 表示对数据进行中心化(均值为0)。
    • pca.scale = F 表示不对数据进行缩放(标准化)。
    • var.shape = 'ColdTolerant' 和 var.col = 'Treatment' 用不同的形状和颜色表示不同的冷耐心和处置惩罚组。
    • var.lab = 'Name' 显示样本的标签(样本名称)。
    • pal.col = 'aaas' 利用指定的颜色调色板。
    • 最后,利用ggsave()生存PCA图。

  • t-SNE分析
    1. p3 = plot_tsne(tm, th, pct.exp = .8, perp = 3, iter = 1200, seed = 2,
    2.     var.shape = 'ColdTolerant', var.col = 'Treatment', var.lab = 'Name',
    3.     legend.pos = 'top.right', legend.dir = 'v', pal.col = 'aaas')
    4. ggsave(file.path(dirw, '11.tsne.pdf'), p3, width = 6, height = 6)
    复制代码

    • plot_tsne()用于执行t-SNE降维,并将高维数据映射到二维空间。
    • pct.exp = .8 利用80%的基因表达数据。
    • perp = 3 设置t-SNE的“狐疑度”参数,这影响数据点的聚集程度。
    • iter = 1200 设置t-SNE迭代的次数。
    • seed = 2 设置随机种子,确保结果可重复。
    • var.shape = 'ColdTolerant' 和 var.col = 'Treatment' 根据样本的冷耐心和处置惩罚组着色与形状区分样本。
    • legend.pos = 'top.right' 设置图例的位置。
    • 最后,生存t-SNE图像。

总结:

这一部门的代码实现了:


  • 层次聚类分析(利用皮尔逊和斯皮尔曼干系系数)。
  • PCA降维分析,可用于观察样本在主成分空间中的分布。
  • t-SNE降维分析,通过非线性方法将高维数据嵌入到二维空间中,帮助明白样本之间的相似性。
这些可视化有助于研究者了解数据的整体结构、样本的聚类情况以及不同处置惩罚组或冷耐心状态下的差别。
好的,接下来是第四部门:数据处置惩罚与生存(过滤样本信息与生存元数据)
第四部门:数据处置惩罚与生存(过滤样本信息与生存元数据)

  1. #{{{ raw: filter/fix samples
  2. th2 = res$th
  3. th2 = complete_sample_list(th2)
  4. fh = file.path(dirw, '01.meta.tsv')
  5. write_tsv(th2, fh, na='')
  6. #}}}
复制代码
解析

  • th2 = res$th
    将之前加载的样本信息 res$th 赋值给 th2。th2 存储的是包含全部样本元数据的数据框。
  • th2 = complete_sample_list(th2)
    complete_sample_list() 是一个自定义的函数,用于过滤或修复样本信息,通常目标是确保每个样本的元数据完整且没有缺失值。这个步骤可以用来清算样本数据,去除不完整或不符合标准的样本。
  • fh = file.path(dirw, '01.meta.tsv')
    创建一个文件路径 fh,用于存储样本的元数据(th2)。路径是在先前定义的工作目次 dirw 下,文件名为 01.meta.tsv。
  • write_tsv(th2, fh, na='')
    将处置惩罚后的 th2 数据框(样本元数据)生存为一个 TSV 文件。write_tsv() 函数将数据框写入文件,并且将缺失值(NA)生存为空字符串。
总结:

这一部门的重要任务是:


  • 过滤样本数据:通过 complete_sample_list() 函数清算样本元数据,确保数据的完整性。
  • 生存元数据:将清算后的样本信息生存为 TSV 格式的文件(01.meta.tsv),以便后续利用或分析。
好的,最后是第五部门:读取RNA-seq数据并重新进行聚类与可视化
第五部门:读取RNA-seq数据并重新进行聚类与可视化

  1. #{{{ read in
  2. res = rnaseq_cpm(yid)
  3. th = res$th; tm = res$tm; tl = res$tl; th_m = res$th_m; tm_m = res$tm_m
  4. th2 = read_xlsx(glue("~/projects/cold/data/samples.xlsx")) %>%
  5.     replace_na(list(Name='?', ColdTolerant='?'))
  6. th = res$th %>% left_join(th2, by="Genotype") %>%
  7.     mutate(lab = str_c(Genotype, Name, Treatment, ColdTolerant, sep='_'))
  8. tm = res$tm %>% filter(SampleID %in% th$SampleID) %>%
  9.     mutate(value=asinh(CPM))
  10. #}}}
  11. #{{{ hclust & tSNE
  12. p1 = plot_hclust(tm, th, pct.exp = .7, cor.opt = 'pearson', var.col = 'Treatment', expand.x = .4)
  13. ggsave(file.path(dirw, '21.hclust.p.pdf'), p1, width = 6, height = 8)
  14. p1 = plot_hclust(tm, th, pct.exp = .7, cor.opt = 'spearman', var.col = 'Treatment', expand.x = .4)
  15. ggsave(file.path(dirw, '21.hclust.s.pdf'), p1, width = 6, height = 8)
  16. p2 = plot_pca(tm, th, pct.exp = .7, pca.center = T, pca.scale = F,
  17.     var.shape = 'ColdTolerant', var.col = 'Treatment', var.lab = 'Name',
  18.     legend.pos = 'bottom.left', legend.dir = 'v', pal.col = 'aaas')
  19. ggsave(file.path(dirw, '21.pca.pdf'), p2, width = 6, height = 6)
  20. p3 = plot_tsne(tm, th, pct.exp = .8, perp = 3, iter = 1200, seed = 2,
  21.     var.shape = 'ColdTolerant', var.col = 'Treatment', var.lab = 'Name',
  22.     legend.pos = 'top.right', legend.dir = 'v', pal.col = 'aaas')
  23. ggsave(file.path(dirw, '21.tsne.pdf'), p3, width = 6, height = 6)
  24. #}}}
复制代码
解析

  • 重新加载RNA-seq数据
    1. res = rnaseq_cpm(yid)
    2. th = res$th; tm = res$tm; tl = res$tl; th_m = res$th_m; tm_m = res$tm_m
    复制代码
    这一行代码再次加载了RNA-seq数据,执行的是 rnaseq_cpm() 函数,通常该函数会返回一系列数据框:

    • th:样本的元数据(例如样本ID、处置惩罚信息等)。
    • tm:基因的表达矩阵(每个基因在各个样本中的表达量)。
    • tl:基因注释信息。
    • th_m 和 tm_m:大概是颠末一些处置惩罚或调整的数据。

  • 重新加载样本信息并归并
    1. th2 = read_xlsx(glue("~/projects/cold/data/samples.xlsx")) %>%
    2.     replace_na(list(Name='?', ColdTolerant='?'))
    3. th = res$th %>% left_join(th2, by="Genotype") %>%
    4.     mutate(lab = str_c(Genotype, Name, Treatment, ColdTolerant, sep='_'))
    5. tm = res$tm %>% filter(SampleID %in% th$SampleID) %>%
    6.     mutate(value=asinh(CPM))
    复制代码

    • 从Excel文件重新读取样本信息,确保样本的“Name”和“ColdTolerant”列中的缺失值被更换为?。
    • 利用left_join()将样本信息th2与res$th按“Genotype”列归并,创建新的lab列,拼接了样本的基因型、名称、处置惩罚信息和冷耐心信息。
    • 筛选出与th中样本ID匹配的基因表达数据tm,并对CPM值进行反双曲正弦变更。

  • 重新进行聚类分析与可视化
    这一部门的分析和可视化步骤与前面相似,都是对基因表达数据进行聚类分析和降维可视化:

    • 层次聚类(hclust):利用皮尔逊和斯皮尔曼干系系数进行层次聚类,并生存为hclust.p.pdf和hclust.s.pdf。
    • PCA分析:通过plot_pca()进行主成分分析,生存为pca.pdf。
    • t-SNE分析:通过plot_tsne()进行t-SNE降维,生存为tsne.pdf。
    这些分析和可视化结果用于展示数据的高维结构,帮助研究人员明白不同样本之间的关系和差别。

总结:

这一部门的代码做了以下几件事情:


  • 重新加载RNA-seq数据,并更新样本的元数据。
  • 对基因表达数据进行处置惩罚(如变更和归并)。
  • 重新进行聚类分析和降维分析,并生存可视化结果。
这与前面处置惩罚的流程雷同,只是数据加载和分析的操作略有不同。重要是为了确保数据的同等性和正确性,随后生成新的聚类与降维结果。

免责声明:如果侵犯了您的权益,请联系站长,我们会及时删除侵权内容,谢谢合作!更多信息从访问主页:qidao123.com:ToB企服之家,中国第一个企服评测及商务社交产业平台。
回复

使用道具 举报

0 个回复

正序浏览

快速回复

您需要登录后才可以回帖 登录 or 立即注册

本版积分规则

海哥

金牌会员
这个人很懒什么都没写!
快速回复 返回顶部 返回列表