ERROR with step of MAPS_PEAK:MAPS_RAW2BG2

[realzhang]

Description of the bug

Dear HiCar authors:
When I run hicar with pilot seq data, i.e., sequenced at low-depth, it errs with the following info. I check the file such as reg_raw.chr1.adultCM.5k.and.MAPS2_pospoisson, and find indeed no "ratio2" column exists. Actually, there are only 2~6 lines in the files naming .5k.and.. ( more than 900 lines in .5k.xor., and "ratio2" column does exist")
I've no idea of the difference between "and and "xor" pospoisson files, and maybe, the low seq depth caused this error?
thank you very much.

create single files from reg_raw to bedgraph for cool load and juicerbox [35/7919]

#########################################
versions <- c("NFCORE_HICAR:HICAR:MAPS_PEAK:MAPS_RAW2BG2:", " MAPS: 1.1.0")
writeLines(versions, "versions.yml") # write versions.yml

options("scipen"=999)
library(data.table)
RESOLUTION = as.numeric(5000)
BG2_OUT = "adultCM_5000.bg2"
GI_OUT = "adultCM_5000.ginteractions"
infs = strsplit("reg_raw.chr1.adultCM.5k.and.MAPS2_pospoisson reg_raw.chr1.adultCM.5k.xor.MAPS2_pospoisson reg_raw.chr10.adultCM.5k.and.MAPS2_pospoisson reg_raw.chr10.adultCM.5k.xor.MAPS2_pospoisson reg_raw.chr11.adultCM.5k.and.MAPS2_pospoisson reg_raw.
chr11.adultCM.5k.xor.MAPS2_pospoisson reg_raw.chr12.adultCM.5k.and.MAPS2_pospoisson reg_raw.chr12.adultCM.5k.xor.MAPS2_pospoisson reg_raw.chr13.adultCM.5k.and.MAPS2_pospoisson reg_raw.chr13.adultCM.5k.xor.MAPS2_pospoisson reg_raw.chr14.adultCM.5k.and.MAPS
2_pospoisson reg_raw.chr14.adultCM.5k.xor.MAPS2_pospoisson reg_raw.chr15.adultCM.5k.and.MAPS2_pospoisson reg_raw.chr15.adultCM.5k.xor.MAPS2_pospoisson reg_raw.chr16.adultCM.5k.and.MAPS2_pospoisson reg_raw.chr16.adultCM.5k.xor.MAPS2_pospoisson reg_raw.chr1
7.adultCM.5k.and.MAPS2_pospoisson reg_raw.chr17.adultCM.5k.xor.MAPS2_pospoisson reg_raw.chr18.adultCM.5k.and.MAPS2_pospoisson reg_raw.chr18.adultCM.5k.xor.MAPS2_pospoisson reg_raw.chr19.adultCM.5k.and.MAPS2_pospoisson reg_raw.chr19.adultCM.5k.xor.MAPS2_po
spoisson reg_raw.chr2.adultCM.5k.and.MAPS2_pospoisson reg_raw.chr2.adultCM.5k.xor.MAPS2_pospoisson reg_raw.chr3.adultCM.5k.and.MAPS2_pospoisson reg_raw.chr3.adultCM.5k.xor.MAPS2_pospoisson reg_raw.chr4.adultCM.5k.and.MAPS2_pospoisson reg_raw.chr4.adultCM.
5k.xor.MAPS2_pospoisson reg_raw.chr5.adultCM.5k.and.MAPS2_pospoisson reg_raw.chr5.adultCM.5k.xor.MAPS2_pospoisson reg_raw.chr6.adultCM.5k.and.MAPS2_pospoisson reg_raw.chr6.adultCM.5k.xor.MAPS2_pospoisson reg_raw.chr7.adultCM.5k.and.MAPS2_pospoisson reg_ra
w.chr7.adultCM.5k.xor.MAPS2_pospoisson reg_raw.chr8.adultCM.5k.and.MAPS2_pospoisson reg_raw.chr8.adultCM.5k.xor.MAPS2_pospoisson reg_raw.chr9.adultCM.5k.and.MAPS2_pospoisson reg_raw.chr9.adultCM.5k.xor.MAPS2_pospoisson", "\s+")[[1]]

peaks_final_out <- lapply(infs, function(inf){
peaks = as.data.table(read.table(inf, header=TRUE, stringsAsFactors=FALSE))
if(nrow(peaks)==0){
peaks$bin1_end <- peaks$bin2_end <- numeric(0)
}else{
peaks$bin1_end = peaks$bin1_mid + RESOLUTION
peaks$bin2_end = peaks$bin2_mid + RESOLUTION
}
peaks_final = subset(peaks, select = c( "chr", "bin1_mid", "bin1_end",
"chr", "bin2_mid", "bin2_end",
"ratio2"))
colnames(peaks_final) = c('chrom1', 'start1', 'end1', 'chrom2', 'start2', 'end2', 'count')
peaks_final
})
peaks_final_out <- do.call(rbind, peaks_final_out)
peaks_final_out <- peaks_final_out[peaks_final_out$count>0, , drop=FALSE]
peaks_final_out <- peaks_final_out[order(peaks_final_out$chrom1,
peaks_final_out$start1,
peaks_final_out$chrom2,
peaks_final_out$start2), , drop=FALSE]
write.table(peaks_final_out, BG2_OUT, row.names = FALSE, col.names = FALSE, quote=FALSE, sep=' ')
#print
peaks_final_out = cbind(0, peaks_final_out[, c('chrom1', 'start1')], 0,
0, peaks_final_out[, c('chrom2', 'start2')], 1,
peaks_final_out[, 'count'])
write.table(peaks_final_out, GI_OUT, row.names = FALSE, col.names = FALSE, quote=FALSE, sep=' ')

Command exit status:
1

Command output:
(empty)

Command error:
Error in [.data.table(x, r, vars, with = FALSE) :
column(s) not found: ratio2
Calls: lapply ... FUN -> subset -> subset.data.table -> [ -> [.data.table
Execution halted

Command used and terminal output

No response

Relevant files

No response

System information

No response

[YCMaCY]

I also encountered this problem before, I found that the number of MAPS “and peak” was too low, when I increased the number of HICAR sequencing reads, it was solved。

[jianhong]

@realzhang Sorry for the late reply. Did you tried to lower the cutoff values for MAPS?