Signac CallPeaks from multiple fragment files
0
3
Entering edit mode
3.3 years ago
LacquerHed ▴ 30

I am attempting to run Macs2 CallPeaks on some multiome data and running into a problem when attempting to run CallPeaks command on multiple fragment file paths in Seurat object.

peaks<-CallPeaks(DataCombined, macs2.path = "/anaconda3/bin/macs2")

FileNotFoundError: [Errno 2] No such file or directory: '/Users/Desktop/multiome/sc291/atac_fragments.tsv.gz /Users/Desktop/multiome/sc292/atac_fragments.tsv.gz /Users/Desktop/multiome/sc293/atac_fragments.tsv.gz'

So basically its listing all three paths together as a single path, not sure how to resolve this. TIA

Seurat Macs2 Signac • 3.8k views
ADD COMMENT
0
Entering edit mode

Hi, I'm having the very same issue, have you solved it?

Thanks! Harris

ADD REPLY
1
Entering edit mode

What is the class of your input data? You might want to include all code thus far.

ADD REPLY
1
Entering edit mode

My input data is a Seurat object, created similarly to here: https://satijalab.org/signac/articles/pbmc_multiomic.html. I merged two such objects, each with a different fragment file:

p0f.inhib <- merge( x = p0f.inhib, y = list(p0f2.inhib), add.cell.ids = c("p0f", "p0f2") )

Of note, my CallPeaks function works perfectly well on each individual object, but not on the merged one. Here is my output:

> peaks <- CallPeaks(p0f.inhib, macs2.path = "/Users/kaplan/MyPythonEnv/bin/MACS2")
INFO  @ Mon, 16 Aug 2021 17:14:22: 
# Command line: callpeak -t /Users/kaplan/Dropbox/Dulac_Lab/Data/Sequencing_Analysis/10x_multiome/july_2021_run2_P0_P10_P18_P28/P0F/atac_fragments.tsv.gz /Users/kaplan/Dropbox/Dulac_Lab/Data/Sequencing_Analysis/10x_multiome/july_2021_run2_P0_P10_P18_P28/P0F2/atac_fragments.tsv.gz -g 2.7e+09 -f BED --nomodel --extsize 200 --shift -100 -n SeuratProject --outdir /var/folders/z6/r0r4szdd68zcjn79vcd0hk_m001q1j/T//RtmpKrGSsm
# ARGUMENTS LIST:
# name = SeuratProject
# format = BED
# ChIP-seq file = ['/Users/kaplan/Dropbox/Dulac_Lab/Data/Sequencing_Analysis/10x_multiome/july_2021_run2_P0_P10_P18_P28/P0F/atac_fragments.tsv.gz /Users/kaplan/Dropbox/Dulac_Lab/Data/Sequencing_Analysis/10x_multiome/july_2021_run2_P0_P10_P18_P28/P0F2/atac_fragments.tsv.gz']
# control file = None
# effective genome size = 2.70e+09
# band width = 300
# model fold = [5, 50]
# qvalue cutoff = 5.00e-02
# The maximum gap between significant sites is assigned as the read length/tag size.
# The minimum length of peaks is assigned as the predicted fragment length "d".
# Larger dataset will be scaled towards smaller dataset.
# Range for calculating regional lambda is: 10000 bps
# Broad region calling is off
# Paired-End mode is off

INFO  @ Mon, 16 Aug 2021 17:14:22: #1 read tag files... 
INFO  @ Mon, 16 Aug 2021 17:14:22: #1 read treatment tags... 
Traceback (most recent call last):
  File "/Users/kaplan/MyPythonEnv/bin/MACS2", line 653, in <module>
    main()
  File "/Users/kaplan/MyPythonEnv/bin/MACS2", line 51, in main
    run( args )
  File "/Users/kaplan/MyPythonEnv/lib/python3.9/site-packages/MACS2/callpeak_cmd.py", line 65, in run
    else:       (treat, control) = load_tag_files_options  (options)
  File "/Users/kaplan/MyPythonEnv/lib/python3.9/site-packages/MACS2/callpeak_cmd.py", line 387, in load_tag_files_options
    tp = options.parser(options.tfile[0], buffer_size=options.buffer_size)
  File "MACS2/IO/Parser.pyx", line 330, in MACS2.IO.Parser.GenericParser.__init__
  File "/Library/Frameworks/Python.framework/Versions/3.9/lib/python3.9/gzip.py", line 58, in open
    binary_file = GzipFile(filename, gz_mode, compresslevel)
  File "/Library/Frameworks/Python.framework/Versions/3.9/lib/python3.9/gzip.py", line 173, in __init__
    fileobj = self.myfileobj = builtins.open(filename, mode or 'rb')
FileNotFoundError: [Errno 2] No such file or directory: '/Users/kaplan/Dropbox/Dulac_Lab/Data/Sequencing_Analysis/10x_multiome/july_2021_run2_P0_P10_P18_P28/P0F/atac_fragments.tsv.gz /Users/kaplan/Dropbox/Dulac_Lab/Data/Sequencing_Analysis/10x_multiome/july_2021_run2_P0_P10_P18_P28/P0F2/atac_fragments.tsv.gz'
Error in file(file, "rt") : cannot open the connection
In addition: Warning message:
In file(file, "rt") :
  cannot open file '/var/folders/z6/r0r4szdd68zcjn79vcd0hk_m001q1j/T//RtmpKrGSsm/SeuratProject_peaks.narrowPeak': No such file or directory
ADD REPLY
0
Entering edit mode

Have not had luck with this yet, you?

ADD REPLY
0
Entering edit mode

Unfortunately not. Will let you know if I find a solution.

ADD REPLY
0
Entering edit mode

FYI - I decided to simply proceed with calling peaks using each fragment file separately.

ADD REPLY
0
Entering edit mode

They added a fix for this in the newest development version of Signac. Works well for me now.

ADD REPLY

Login before adding your answer.

Traffic: 2586 users visited in the last hour
Help About
FAQ
Access RSS
API
Stats

Use of this site constitutes acceptance of our User Agreement and Privacy Policy.

Powered by the version 2.3.6