converting HiC output files to juicebox input
Shane
nservant
Indeed, it looks like the sort does not work for memory issue. But this is a space issue in your root tmp, ie. /tmp/, not ./tmp/
Could you try ;
awk '{$4=$4!="+"; $7=$7!="+"; split($9, frag1, "_"); split($10, frag2, "_"); } $2<=$5{print $1, $4, $2, $3, frag1[3], $7, $5, $6, frag2[3], $11, $12 }$5<$2{ print $1, $7, $5, $6, frag2[3], $4, $2, $3, frag1[3], $12, $11}' hic_output/data/input_1/input_1_allValidPairs > ./tmp/test.out
sort -k3,3d -k7,7d -S 90 ./tmp/test.out > ./tmp/test_allValidPairs.pre_juicebox_sorted
These are the commands run by hicpro2juicebox.
I think that the sort will crash ...
N
Shane
nservant
Does it make sense ?
Shane
Now the task is still running, I will post the results here once it is finished.
Shane
nservant
I think that you are using Juicebox ? and not the Juicebox command line tools (juicebox_clt) ?
You need to use the command line tools to generate the .hic file.
N
Shane
Mandy Wong
Generating Juicebox input files ...
nservant
Yes it can be long. From my experience, there are two things you can try to do.
1/ remove the -r option. It is useful if you want to see your maps at the restriction fragment resolution, but it also requires a huge sequencing depth.
2/ could check the list of chromosome in the /chrom_hg19.sizes file ? I realized that the more chromosomes you have, the longer it is.
I fyou have any 'random' chromosomes, I would suggest to remove them.
N
Mandy Wong
awk '{$4=$4!="+"; $7=$7!="+"; n1=split($9, frag1, "_"); n2=split($10, frag2, "_"); } $2<=$5{print $1, $4, $2, $3, frag1[n1],
$7, $5, $6, frag2[n2], $11, $12 }$5<$2{ print $1, $7, $5, $6, frag2[n2], $4, $2, $3, frag1[n1], $12, $11}' $VALIDPAIRS | sort -
k3,3d -k7,7d -S 90 > ${TEMP}/$$_allValidPairs.pre_juicebox_sorted
and the .prejuicebox_sorted file is still empty.(0K)It seems like making changes (1) & (2) you've suggested above still yields the same results. Is there something wrong maybe regarding the format of the output for allValidPairs that could be causing the problem?
I attach here the head of my .allValidPairs file:
SRR1658673.259994198 chr1 10776 + chr7 131533799 + 570 HIC_chr1_1 HIC_chr7_319776 0 42
SRR1658673.72579757 chr1 12259 + chr15 93563685 - 341 HIC_chr1_2 HIC_chr15_188640 2 40
SRR1658676.125252832 chr1 13028 + chr4 134982856 + 426 HIC_chr1_6 HIC_chr4_308064 30 42
SRR1658676.120863062 chr1 13028 - chr15 102348304 - 454 HIC_chr1_6 HIC_chr15_209860 31 42
SRR1658676.152863214 chr1 13028 - chr18 77468066 - 361 HIC_chr1_6 HIC_chr18_177304 30 42
SRR1658673.90423160 chr1 13028 - chr19 560618 + 439 HIC_chr1_6 HIC_chr19_1329 0 2
SRR1658673.201320596 chr1 13028 - chr19 16679875 + 297 HIC_chr1_6 HIC_chr19_54907 31 42
SRR1658676.103541592 chr1 13029 + chr3 78354693 - 456 HIC_chr1_6 HIC_chr3_194482 6 23
SRR1658676.219111220 chr1 13029 - chr11 111314436 - 393 HIC_chr1_6 HIC_chr11_265725 31 42
SRR1658676.81193276 chr1 13031 + chr15 86012106 + 365 HIC_chr1_6 HIC_chr15_168511 11 42
and the file is 70GB....
nservant
Yes indeed. They move the command line tools into Juicer now.
So the jar should be this one ;
https://github.com/theaidenlab/juicer/wiki/Juicer-Tools-Quick-Start
Regarding you issu, I'm a bit confused by the fact that the awk command did not finish.
My feeling is that this is more related to the sort of a 70Gb file.
To double check that, could you please run the same command on the head of your file ?
And then the same on your entire file, but without the sort ?
Thanks
Mandy Wong
Problem with creating fragment-delimited maps, NullPointerException.
This could be due to a null fragment map or to a mismatch in the chromosome name in the fragment map vis-a-vis the input file or chrom.sizes file.
Exiting.
done !
Does that mean some possible formatting error in my file?
nservant
It means that it passes through the awk command.
Did you try to run only the awk command on your valid pairs file, without the sort.
My feeling is that the sort of the 70Gb file is the issue ...
Mandy Wong
nservant
Otherwise, how many RAM do you have ? Could you have access to a more powerfull machine ?
Mandy Wong
nservant
One idea would be to specify the TMP_DIR of sort using -T option.
I think that by default, sort writes its temporary files in your current folder. If for any reason you do not have enough disk, it could crash.
Let me know if I can help.
If necessary, you can send me your file in a private message. I can try here.
Best wishes
Mandy Wong
${TEMP}/$$_allValidPairs.pre_juicebox_sorted <- which is be in ./tmp , where . is whatever the directory I ran the script. Is it odd that, when I include the sort part of the awk command, the script doesn't show signs of crashing(quitting and going back to the command prompt), and yet the pre_juicebox_sorted file is 0KB the whole time? or is it that if the sorting was never finished and the program has crashed, this is what we expect to see?
Thank you so much for you help on running the script and the whole pipeline :)
nservant
And indeed, my feeling is that the sort never writes its output.
So finally, did it run without the sort ?
Mandy Wong
nservant
Could you try to remove the -S 90 ? using this option, it allows sort to use 90% of your RAM ... it might be too much.
N
Mandy Wong
awk '{$4=$4!="+"; $7=$7!="+"; n1=split($9, frag1, "_"); n2=split($10, frag2, "_"); } $2<=$5{print $1, $4, $2, $3, frag1[n1], $7, $5, $6, frag2[n2], $11, $12 }$5<$2{ print $1, $7, $5, $6, f rag2[n2], $4, $2, $3, frag1[n1], $12, $11}' $VALIDPAIRS > ./tmp/test.out
##sort -k3,3d -k7,7d -S 90 ./tmp/test.out > ${TEMP}/$$_allValidPairs.pre_juicebox_sorted
##awk '{$4=$4!="+"; $7=$7!="+"; n1=split($9, frag1, "_"); n2=split($10, frag2, "_"); } $2<=$5{print $1, $4, $2, $3, frag1[n1], $7, $5, $6, frag2[n2], $11, $12 }$5<$2{ print $1, $7, $5, $6, f rag2[n2], $4, $2, $3, frag1[n1], $12, $11}' $VALIDPAIRS > ./tmp/test.out
##sort -k3,3d -k7,7d ./tmp/test.out > ${TEMP}/$$_allValidPairs.pre_juicebox_sorted
pengt...@gmail.com
在 2018年4月6日星期五 UTC+8上午3:02:44，nservant写道：
nservant
N
