我有一個 R 資料框,其中包含 SNP 對(染色體編號和堿基對位置在不同的列中)和每對的 p 值。資料如下所示:
chr_snp1 loc_snp1 chr_snp2 loc_snp2 pval
1 278474050 2 57386477 7.43e-08
1 275620856 2 57386662 1.08e-07
1 144075771 3 109909704 1.02e-06
1 144075771 3 111344453 2.06e-06
2 103701229 7 56369738 3.83e-08
2 102990566 7 56407691 1.07e-07
我想成對洗掉冗余。loc_snp1但是,冗余被定義為與由和loc_snp2列以及chr_snp1和中的值確定的另一個相似的對chr_snp2。loc_snp1我想對和中的值施加的數值限制loc_snp2是 10,000,000,對于同一染色體組合上的對。該資料按染色體(chr_snp1和chr_snp2)排序,然后是堿基對位置(loc_snp1和loc_snp2),然后按p值(pval)排序。
基本上我希望腳本檢查 in的值是否在上一行loc_snp1中的值的 10,000,000 范圍內,以及它們是否在. 然后,還要檢查 in的值是否在上一行 in 的值的 10,000,000 以內,以及它們是否具有相同的 in 值。如果所有四個條件都為真,則只保留 pval 最低的行并丟棄另一個。loc_snp1chr_snp1loc_snp2loc_snp2chr_snp2
換句話說,只保留最好的一對并洗掉所有其他接近的(具有相同的染色體組合并且在它們各自的染色體上的其他對的 10,000,000 個堿基對內)。
這將導致資料框看起來像:
chr_snp1 loc_snp1 chr_snp2 loc_snp2 pval
1 278474050 2 57386477 7.43e-08
1 144075771 3 109909704 1.02e-06
2 103701229 7 56369738 3.83e-08
當然,腳本不必實際檢查上面的行。我認為有更優雅的方法來解決這個問題。
uj5u.com熱心網友回復:
這是選項。如果您已經知道 chr_snp1 和 chr_snp2 需要相同,則可以先將它們分組,然后計算差異。誠然,如果您在資料集中的其他地方有重復的 chr_snp1,這可能會失敗。
library(tidyverse)
df |>
group_by(chr_snp1, chr_snp2) |>
mutate(grp = (abs(loc_snp1 - lag(loc_snp1, default = first(loc_snp1))) < 10e6) &
(abs(loc_snp2 - lag(loc_snp2, default = first(loc_snp2))) < 10e6)) |>
group_by(grp, .add=TRUE) |>
filter(pval == min(pval)) |>
ungroup()|>
select(-grp)
#> # A tibble: 3 x 5
#> chr_snp1 loc_snp1 chr_snp2 loc_snp2 pval
#> <dbl> <dbl> <dbl> <dbl> <dbl>
#> 1 1 278474050 2 57386477 0.0000000743
#> 2 1 144075771 3 109909704 0.00000102
#> 3 2 103701229 7 56369738 0.0000000383
uj5u.com熱心網友回復:
df <- textConnection('
chr_snp1, loc_snp1, chr_snp2, loc_snp2, pval
1, 278474050, 2, 57386477, 7.43e-08
1, 275620856, 2, 57386662, 1.08e-07
1, 194406449, 2, 182907442, 2.22e-06
1, 144075771, 3, 109909704, 1.02e-06
1, 144075771, 3, 111344453, 1.02e-06
2, 103701229, 7, 56369738, 1.02e-06
2, 102990566, 7, 56407691, 1.02e-06
') |> read.csv(header = TRUE)
# filter based on first 4 criteria
threshold <- 1e7
idx <- (
((diff(df$loc_snp1) |> abs()) <= threshold) &
(diff(df$chr_snp1) == 0) &
((diff(df$loc_snp2) |> abs()) <= threshold) &
(diff(df$chr_snp2) == 0)
) |> which()
df <- df[ idx, ]
# now retain only the pair with lowest p-value
df <- split(df, paste0(df$chr_snp1, df$chr_snp2)) |>
lapply(function(x) {
idx <- which.min(x$pval)
x[ idx, ]
}) |>
data.table::rbindlist() |>
as.data.frame()
# now retain only the pair with lowest p-value (mirror redundancies)
df <- split(df, paste0(df$chr_snp2, df$chr_snp1)) |>
lapply(function(x) {
idx <- which.min(x$pval)
x[ idx, ]
}) |>
data.table::rbindlist() |>
as.data.frame()
print(df)
chr_snp1 loc_snp1 chr_snp2 loc_snp2 pval
1 1 278474050 2 57386477 7.43e-08
3 1 144075771 3 109909704 1.02e-06
5 2 103701229 7 56369738 1.02e-06
uj5u.com熱心網友回復:
在 AndS 的初始修剪步驟之后,一位朋友能夠撰寫如何處理鏡像重復的代碼。
Cart_Inter_Pruned$skip = FALSE
df = NULL
threshold = 10e7
for(i in 1:nrow(Cart_Inter_Pruned)){
if(Cart_Inter_Pruned[i,]$skip) next
ordered_key1 = sort(array(unlist(Cart_Inter_Pruned[i,c(1,3)])))
for(j in (i 1):nrow(Cart_Inter_Pruned)){
if(Cart_Inter_Pruned[j,]$skip) next
ordered_key2 = sort(array(unlist(Cart_Inter_Pruned[j,c(1,3)])))
if(paste(ordered_key1,collapse=",") == paste(ordered_key2,collapse = ",")){
if_switch = ordered_key1[1] == ordered_key2[2]
print(paste(i,j,ordered_key1))
if(if_switch){
if(abs(Cart_Inter_Pruned[i,2]- Cart_Inter_Pruned[j,4]) <= threshold && abs(Cart_Inter_Pruned[i,4] - Cart_Inter_Pruned[j,2]) <= threshold){
Cart_Inter_Pruned[j,]$skip = TRUE
Cart_Inter_Pruned[i,5] = min(as.numeric(Cart_Inter_Pruned[i,5] ,as.numeric(Cart_Inter_Pruned[j,5])))
}
}else{
if(abs(Cart_Inter_Pruned[i,2]- Cart_Inter_Pruned[j,2]) <= threshold && abs(Cart_Inter_Pruned[i,4] - Cart_Inter_Pruned[j,4]) <= threshold){
Cart_Inter_Pruned[j,]$skip = TRUE
Cart_Inter_Pruned[i,5] = min(as.numeric(Cart_Inter_Pruned[i,5] ,as.numeric(Cart_Inter_Pruned[j,5])))
}
}
}
}
}
轉載請註明出處,本文鏈接:https://www.uj5u.com/qukuanlian/526994.html
標籤:r数据框条件语句数据操作
