xryanglab / RiboMiner

a toolset for mining multi-dimensional features of the translatome with ribosome profiling data
https://pubmed.ncbi.nlm.nih.gov/32738892/
GNU General Public License v3.0
13 stars 7 forks source link

Enrichment analysis for single transcript error #6

Closed junjunlab closed 2 years ago

junjunlab commented 2 years ago

Dear author, thank you very much for this great software to analysis Ribo-seq data. Recently I use this software to analysis seRP Ribo-seq enrichment, I want to make a single trancript enrichment ,this is my code:

EnrichmentAnalysis --ctrl 2.RiboDensity/density_bc_cds_codon_density.txt --treat 2.RiboDensity/density_IgG_cds_codon_density.txt -c longest.transcripts.info.txt -o 4.single-enrich/IgG_M1 -U codon -M RPKM -l 150 -n 10 -m 1 -e 30 --CI 0.95 -u 0 -d 500 -S M1.txt --id-type transcript_id

but I got a error like this:

19693 transcripts will be used in the follow analysis.

There are 1 transcripts from M1.txt used for following analysis. There are 1 transcripts with both start and stop codon will be used for following analysis. Traceback (most recent call last): File "/home/zhoulab/anaconda3/envs/ribominer/bin/EnrichmentAnalysis", line 10, in sys.exit(main()) File "/home/zhoulab/anaconda3/envs/ribominer/lib/python3.7/site-packages/RiboMiner/EnrichmentAnalysis.py", line 295, in main min_cds_codon,min_cds_counts,downstream_codon,upstream_codon,norm_exclude_codon,min_norm_region_counts,mode,unit,confidence) File "/home/zhoulab/anaconda3/envs/ribominer/lib/python3.7/site-packages/RiboMiner/EnrichmentAnalysis.py", line 150, in enrichment_ratio startratio[terms]=np.mean(startNormedWindowsList[np.where(startWindowsExistPosList[:,terms]==1),terms]) IndexError: too many indices for array: array is 1-dimensional, but 2 were indexed

So I don't know what's the problem, I am looking forward to your reply!

sherkinglee commented 2 years ago

Hi, Thanks for your interest in RiboMiner. Could you tell me which single transcript you want to analyze? I will test it and tell you the solution.

sherkinglee commented 2 years ago

Hi, I have tested it with a single transcript or a subset of transcripts. The problem may be that your single transcript does not meet the criteria with "-l 150 -n 10 -m 1" which means the transcript has to be longer than 150 codons and with RPKM>=10. It may work if you set "-l 0 -n 0 -m 0" which means you do not set any criteria for this transcript. And actually, the EnrichmentAnalysis function was designed for metagene analysis with a subset of transcripts. If you want to study whether a protein could be co-translated with another transcript, I recommend you use the command as follows:

## output all codon ratio values for all transcripts without any filtering
EnrichmentAnalysis --ctrl 2.RiboDensity/density_bc_cds_codon_density.txt --treat 2.RiboDensity/density_IgG_cds_codon_density.txt -c longest.transcripts.info.txt -o 4.single-enrich/all -U codon -M RPKM -l 0 -n 0 -m 0 -e 30 --CI 0.95 -u 0 -d 500

## and if your interested transcript has a relative abundance, it would be left in the output file "4.single-enrich/all_codon_ratio.txt". Then plot the enrichment ratio for this transcript:
EnrichmentAnalysisForSingleTrans -i 4.single-enrich/all_codon_ratio.txt -s GAPDH -o GAPDH -c longest.trans.info.txt  --id-type gene_name  --slide-window y --axhline 1

I hope this could help you solve this problem. If you still have any questions, please feel free to let me know.

Best, Sherking

junjunlab commented 2 years ago

Thank you very much for your reply! I will have a try!