ay-lab / fithic

Fit-Hi-C is a tool for assigning statistical confidence estimates to chromosomal contact maps produced by genome-wide genome architecture assays such as Hi-C.
MIT License
79 stars 16 forks source link

Running error #32

Closed JhinAir closed 4 years ago

JhinAir commented 4 years ago

Hi, I got a such error when running fithic with command: python /home/user/liu/Software/fithic/fithic/fithic.py -f fithic.fragmentMappability.gz -i fithic.interactionCounts.gz -o ./Inter -r 40000 -t fithic.biases.gz -x interOnly Btw, the test running is working successfully.

####################################################### Reading fragments file from: fithic.fragmentMappability.gz Reading interactions file from: fithic.interactionCounts.gz Output path created ./Inter Fixed size option detected... Fast version of FitHiC will be used Resolution is 40.0 kb Reading bias file from: fithic.biases.gz The number of spline passes is 1 The number of bins is 100 The number of reads required to consider an interaction is 1 The name of the library for outputted files will be Emu Upper Distance threshold is inf Lower Distance threshold is 0 Only inter-chromosomal regions will be analyzed Lower bound of bias values is 0.5 Upper bound of bias values is 2 All arguments processed. Running FitHiC now...

Reading the contact counts file to generate bins... Interactions file read. Time took 154.33874917030334 Traceback (most recent call last): File "/home/user/liu/Software/fithic/fithic/fithic.py", line 1317, in main() File "/home/user/liu/Software/fithic/fithic/fithic.py", line 323, in main (binStats,noOfFrags, maxPossibleGenomicDist, possibleIntraInRangeCount, possibleInterAllCount, interChrProb, baselineIntraChrProb) = generate_FragPairs(observedInterAllCount, observedInterAllSum, binStats, fragsFile, resolution) File "/home/user/liu/Software/fithic/fithic/fithic.py", line 597, in generate_FragPairs maxFrags[ch]=max([int(i)-resolution/2 for i in allFragsDic[ch]]) ValueError: max() arg is an empty sequence #######################################################

Could you please help with it?

Best~ Jing

ay-lab commented 4 years ago

Hi, seems like there may be a problem with the format of your input files. The dictionary structure supposed to have a list of fragments in it per each chromosome seems to be empty. This may be a problem with your fragment mappability file

JhinAir commented 4 years ago

Hi, thanks for the reply. The fragment file seems to be ok, as following: ######################## 1 0 20000 621 1 1 40000 60000 1073 1 1 80000 100000 837 1 1 120000 140000 976 1 1 160000 180000 1793 1 1 200000 220000 846 1 ######################## The inputs were produced by HiCPro2FitHiC.py.

Best, Jing

ay-lab commented 4 years ago

what I can guess is that you have a chromosome (maybe chr Y if human) in your input file that has no fragments above mappability threshold so hence the error. Try removing chrs that are unmappable. If it doesn't work please share your input files. thanks

JhinAir commented 4 years ago

Hi,

I've checked that all chromosomes have the mappability threshold. The inputs are too big (about 50Mb), so I put them on the google drive. https://drive.google.com/drive/folders/1_wm-CAsXY2oDzRGPLEuPtIIIT5UPIhqA. Please let me know if they are not available. Thanks a lot!

Jing

ay-lab commented 4 years ago

as I have suspected, there are scaffolds in your fragment list and bias file that do not have any valid/mappable fragment. See the examples below. if you don't remove them fithic won't work. Scaffold_274 20000 -1 Scaffold_275 20000 -1 Scaffold_275 60000 -1 Scaffold_275 100000 -1 Scaffold_275 140000 -1 Scaffold_275 180000 -1

JhinAir commented 4 years ago

It's working! Thanks very much! I thought fithic would deal with those -1 value. I have another question that how the ExpCC is calculated for inter-chr interaction, as I found some ExpCC values are equal to zero?

ay-lab commented 4 years ago

It is in the below piece of code. If bias values are outside the desired range then it gets set to zero.

        if allReg or interOnly:
            prior_p=interChrProb*(bias1*bias2)
            p_val=scsp.bdtrc(interxn.getCount()-1,observedInterAllSum,prior_p)
            interCount += 1
            # computing expected contact count
            if ((bias1 >= biasLowerBound) and (bias1 <= biasUpperBound) and (bias2 >= biasLowerBound) and (bias2 <= biasUpperBound)):
                expected_CC = (observedInterAllSum * prior_p)
            else:
                expected_CC = 0