网站建设参考论文,做家教网站代理,saas建站没有网站源代码么,媒易网络网站建设培训变异位点染色体分布图 今天分享的内容是通过RIdeogram包绘制染色体位点分布图#xff0c;并介绍一种展示差异位点的方法。 在遗传学研究中#xff0c;通过测序等方式获得了基因组上某些位置的基因型信息。 如下表#xff0c;第一列是变异位点的ID#xff0c;第二列是染色体… 变异位点染色体分布图 今天分享的内容是通过RIdeogram包绘制染色体位点分布图并介绍一种展示差异位点的方法。 在遗传学研究中通过测序等方式获得了基因组上某些位置的基因型信息。 如下表第一列是变异位点的ID第二列是染色体第三列是物理位置最后两列是两个样品的基因型。 这个文件通常有好几万行要想快速从中获得有效信息可以借助变异位点染色体分布图功能实现。比如下图红色部分表示样品A和样品B差异的位置一眼就能看出哪些位置比较关键。 数据准备 主要的输入文件有两个第一个是染色体结构信息长度、着丝粒位置等第二个是待分析标注的位点信息不同材料在某位置的基因型如下所示 操作方法 首先安装并加载R包 library(RIdeogram)library(tidyverse)library(xlsx) 设置参数 将整理好的文件保存为csv格式输入文件名。centwidth是着丝粒显示宽度singewidth是单个位点的显示宽度这两个值可以自定义设置。如果设置的值越大那么线条就越粗 my_filename - data.csvmy_centwidth - 5000000my_color - c(#e0e0e0, #304ffe, #dd2c00,#00bfa5)my_singewidth - 5000000 添加染色体信息 由于绘图时需要每个染色体的长度起始和终止位置和着丝粒的位置因此提前准备一个txt文档按如下格式整理数据用于后续绘图。 下面是读取染色体信息的代码 df_chr_pos_cent - read.table(./Ref_chromedata.txt,header T)df_chr_pos_cent - cbind(df_chr_pos_cent$chr,0, df_chr_pos_cent$length, df_chr_pos_cent$centromere, as.numeric(df_chr_pos_cent$centromere) my_centwidth)df_chr_pos_cent - as.data.frame(df_chr_pos_cent)df_chr_pos_cent - df_chr_pos_cent[1:21,] %% mutate_at(vars(V2,V3,V4,V5),as.integer)colnames(df_chr_pos_cent) - c(Chr,Start,End,CE_start,CE_end) head(df_chr_pos_cent) Chr Start End CE_start CE_end1 1A 0 594102056 213000000 2180000002 1B 0 689851870 240600000 2456000003 1D 0 495453186 170000000 1750000004 2A 0 780798557 326650000 3316500005 2B 0 801256715 347850000 3528500006 2D 0 651852609 268450000 273450000 先从名为 Ref_chromedata.txt 的文件中读取数据并将其存储在名为 df_chr_pos_cent 的数据框中再选取某些列的内容并将它们拼接成一个新的矩阵自动计算每条染色体的起止位置和着丝粒的起止位置。 创建判定子代和亲本关系的函数 为了更方便的计算某个位点在亲本和子代中的存在情况快速获得子代中突变信息的来源创建一个亲缘关系判定函数。 # 函数输入两个亲本和子代基因型判断子代来源determine_progeny - function(p1,p2,prog,name_p1,name_p2){ if (p1 p2){ if (p1 prog p2 prog){ return(Both) }else{ return(Unknown) } }else{ if (p1 prog){ return(name_p1) }else{ if (p2 prog){ return(name_p2) }else{ return(Unknown) } } }} 这个函数叫做 determine_progeny它用于确定两个亲本p1 和 p2以及一个后代prog之间的关系。 参数说明 p1 - 第一个亲本的标识。 p2 - 第二个亲本的标识。 prog - 后代的标识。 name_p1 - 第一个亲本的名字 name_p2 - 第二个亲本的名字 功能介绍 首先检查如果第一个亲本 p1 和第二个亲本 p2 相等那么它们之间的关系有两种情况 如果它们都等于后代 prog那么返回 Both表示两个亲本都与后代相关。 另外一种情况返回 Unknown表示无法确定它们之间的关系。 如果第一个亲本 p1 和第二个亲本 p2 不相等那么它们之间的关系也有两种情况 如果第一个亲本 p1 等于后代 prog那么返回 name_p1表示第一个亲本是后代的亲本。 如果第二个亲本 p2 等于后代 prog那么返回 name_p2表示第二个亲本是后代的亲本。 最后如果以上条件都不满足返回 Unknown表示无法确定它们之间的关系。 这个函数的目的是根据输入的亲本和后代信息确定它们之间的关系并返回一个相应的标注信息。 创建判定杂合子的函数 如果某些位点为杂合状态需要对其进行识别标注以下提供一种算法进行识别 # 函数判断杂合子decide_hybrid - function(genetype){ tem - str_split(genetype,) if (length(tem[[1]]) ! 2){ return(Error) }else{ first - tem[[1]][1] sencend - tem[[1]][2] if (first sencend){ return(equal) }else{return(diff)} }} 这段代码定义了一个名为decide_hybrid的函数其作用是接受一个字符串参数genetype然后将该字符串拆分成两个字符并检查这两个字符是否相等。 如果输入字符串的长度不等于2函数将返回Error表示输入错误。 否则如果两个字符相等函数将返回equal表示相等。 如果两个字符不相等函数将返回diff表示不相等。这个函数主要用于比较输入的基因型是否相等或不等。 计算每个位点的变异类型 for (i in 1:nrow(df)){ if (df$chr[i] #N/A){ df$type[i] - del next } if (df$pos[i] #N/A){ df$type[i] - del next } if (decide_hybrid(df[i,3]) diff | decide_hybrid(df[i,4]) diff | decide_hybrid(df[i,5]) diff){ df$type[i] - del next } df$type[i] - determine_progeny(df[i,3],df[i,4],df[i,5], colnames(df)[3],colnames(df)[4])} 这段代码的主要作用是根据一系列条件来确定df 数据框中每一行的type列的值并根据不同条件进行相应的处理。如果某些条件满足它将类型设置为 del否则它会计算新的类型并将其存储在数据框中。 绘图数据初步整理 df_marker - df %% filter(type ! del)df_marker$chr - sub(chr,,df_marker$chr)df_marker$Start - df_marker$posdf_marker$End - as.numeric(df_marker$pos) my_singewidthdf_marker$Value[which(df_marker$type Both)] - 4df_marker$Value[which(df_marker$type colnames(df)[3])] - 3df_marker$Value[which(df_marker$type colnames(df)[4])] - 2df_marker$Value[which(df_marker$type Unknown)] - 1 首先从数据框df中筛选出type列不等于del的行然后对筛选后的数据框df_marker进行一系列操作 移除df_marker$chr列中的chr前缀将结果存储回df_marker$chr中。 将df_marker$pos列的值赋给新的列df_marker$Start以便将位置信息重命名。 根据my_singewidth的值计算新的列df_marker$End该列的值等于df_marker$pos加上my_singewidth的数值。 根据条件对df_marker$Value列进行更新将根据type列的不同值分配不同的数值以下为对应规则 如果type等于Both则将df_marker$Value设置为4。如果type等于数据框中的第三列的列名则将df_marker$Value设置为3。如果type等于数据框中的第四列的列名则将df_marker$Value设置为2。如果type等于Unknown则将df_marker$Value设置为1。 主要目的是根据条件对df数据框进行筛选和转换操作创建一个数据框df_marker其中包括重命名列、计算新列以及根据不同的type值赋予相应的数值。 转换数据格式 df_marker_plot - cbind(df_marker[,c(1,7,8,9)])df_marker_plot - df_marker_plot %% mutate_at(vars(Start,End,Value),as.integer)df_marker_plot - df_marker_plot[which(df_marker_plot$chr !Un),]colnames(df_marker_plot) - c(Chr,Start,End,Value) 这段代码的功能是对名为df_marker的数据框进行处理首先选择其中的特定列第1、7、8、9列然后将这些列的数据类型转换为整数类型接着从中删除chr列值不等于Un的行最后重命名这些列为Chr、Start、End和Value用于后续的数据可视化或分析。 绘制染色体全景图 # 绘制全景图ideogram(karyotype df_chr_pos_cent, overlaid df_marker_plot, colorset1 my_color)convertSVG(chromosome.svg, device pdf) 最后使用如上代码即可生成最终的结果图片。该R包非常巧妙的使用SVG语法构建图形并且提供了svg转pdf、png等多种方式。 提示分类型离散数据通过数值映射转换为连续型便可以通过函数进行绘图计算如果有多个维度的数据需要展示可以添加其他类型的标注。 多组学应用 差异表达基因的分布RNA-seq 开放染色质的分布ATAC-seq CTCF结合位点ChIP-seq 突变位点在染色体上的分布WGS DNA甲基化的分布WGBS 外显子在染色体上的分布WES 遗传学应用 变异位点全景图绘制 标记物理位置展示图绘制 遗传多样性来源展示 基因连锁定位结果展示 参考资料https://github.com/cran/RIdeogramhttp://doi.org/10.7717/peerj-cs.251https://www.jianshu.com/p/07ae1fe18071 Tips本文所有示例数据均随机生成不具有任何意义 本文由 mdnice 多平台发布