在excel或R中; 如何在染色体长度上表示1或空白的密度(即SNP是否存在)?

我有一个沿着染色体的物理位置的列,从0.1 – 25,526,585,沿着这条染色体,我有一个1,如果有一个单核苷酸多态性和一个空白,如果没有。 我想制作一个线图或其他东西来显示SNP的任何峰值。我不能手动计算每1000个位置的SNP,因为染色体长度> 25m,位置不规则。 有没有人有一个明智的想法,如何做到这一点,我将非常感激。 数据布局:

Phys_Position Mutant_SNP 0.0 0.1 0.1 0.1 0.1 1 0.1 0.1 1 0.1 0.1 0.2 0.2 0.2 1 0.2 1 0.2 0.2 0.3 0.3 0.7 0.7 0.7 0.7 0.7 0.7 1.4 1.5 1.6 1.7 1.7 1 1.8 1.8 1.9 1.9 2.0 5.4 ... 25,526,585 

输出输出:

 structure(list(PHYS_POS. = c(37, 55, 89, 102, 105, 107, 116, 117, 121, 166), Phys_Position = c(" 0.0 ", " 0.1 ", " 0.1 ", " 0.1 ", " 0.1 ", " 0.1 ", " 0.1 ", " 0.1 ", " 0.1 ", " 0.2 " ), Mutant_SNP = c(NA, NA, NA, NA, 1L, NA, 1L, NA, NA, NA), X = c(NA, NA, NA, NA, NA, NA, NA, NA, NA, NA), X.1 = c(NA, NA, NA, NA, NA, NA, NA, NA, NA, NA), X.2 = c(NA, NA, NA, NA, NA, NA, NA, NA, NA, NA), X.3 = c(NA, NA, NA, NA, NA, NA, NA, NA, NA, NA), X.4 = c(NA, NA, NA, NA, NA, NA, NA, NA, NA, NA), X.5 = c(NA, NA, NA, NA, NA, NA, NA, NA, NA, NA), X.6 = c(NA, NA, NA, NA, NA, NA, NA, NA, NA, NA)), .Names = c("PHYS_POS.", "Phys_Position", "Mutant_SNP", "X", "X.1", "X.2", "X.3", "X.4", "X.5", "X.6"), row.names = c(NA, 10L), class = "data.frame") 

@nyc新数据:

 first.bp mutations 1 0 1000001 0 2000001 0 3000001 0 4000001 0 5000001 0 6000001 0 7000001 0 8000001 0 9000001 0 10000001 0 11000001 0 

下面将返回一个包含每个区间的SNP个数的allSNP 。 您可以将step调整为您想要的任何值,并且您应该根据自己的数据进行操作,而且您还是可以走的。 你要做的最后一件事是做一个结果图。

 #Create Data Phys_Position <- c(0, 5, 10005, 20001) Mutant_SNP <- c(1, 1, 0, 1) df <- data.frame(Phys_Position, Mutant_SNP) df$Phys_Position <- as.numeric(df$Phys_Position) #<-------- added after edit #find first and last value start <- df$Phys_Position[1] limit <- df$Phys_Position[nrow(df)] #intiliaze values step <- 10000 end <- start + step allSNP <- NULL while (start < limit) { subsetData <- subset(df, Phys_Position >= start & Phys_Position < end) nrSNP <- sum(subsetData$Mutant_SNP, na.rm = TRUE) allSNP <- rbind(allSNP, nrSNP) start <- end end <- start + step } 

或者,我们可以使用sapply

 step.size <- 100 pos <- sapply(seq(1, tail(dat$PHYS_POS., 1), step.size), FUN=function(x) sum(dat$Mutant_SNP[dat$PHYS_POS. >= x & dat$PHYS_POS.< (x + step.size - 1)], na.rm = TRUE)) pos [1] 0 2 

根据@ MarcelG的评论,这是对代码的解释。

sapply函数使用sapply从1到最后一个值的序列。 列偏移一个步长。 (对于计算,Phys_Position列中的划分失去了精度,并不是必需的,虽然对于人来说读起来肯定更容易)。然后函数将序列值一次一个地送到在FUN=定义的函数作为variablesx 。 我们使用variablesx的子集,在Mutant_SNP列中的值进行求和。 na.rm = TRUE指定忽略缺失值。

结果可以转换为data.frame并绘制。

 res <- data.frame(first.bp = seq(1, tail(dat$PHYS_POS., 1), step.size), mutations = pos) plot(1, xlim = c(0, max(pos)), ylim = c(0, nrow(res)), type = "n") apply(res, 1, FUN=function(x) segments(0, x[1]/step.size, x[2]))