error using coverageplot() on integrated ATAC data

[afcmalone]

I have been successful integrating 3 ATACseq samples using your vignettes and other posts on github. However, using the final integrated object I get an error when using coverageplot(integrated) This occurs when i try any of the 'assays' in the object but i am assuming the 'mergedpeaks' assay should be used for this.

Here is my code with errors:

DefaultAssay(integrated) <- 'mergedpeaks' CoveragePlot(

• object = integrated,
• region = "chr2-87011729-87035519"
• ) Error in UseMethod(generic = "Fragments", object = object) : no applicable method for 'Fragments' applied to an object of class "Assay"

DefaultAssay(integrated) <- 'peaks' CoveragePlot(

• object = integrated,
• region = "chr2-87011729-87035519"
• ) Error in colnames<-(*tmp*, value = start(x = region):end(x = region)) : attempt to set 'colnames' on an object with less than two dimensions

Thank you in advance.

[timoast]

Looks like mergedpeaks is a standard Assay rather than ChromatinAssay. You can only run CoveragePlot on a ChromatinAssay

[afcmalone]

Thank yo for your quick reply.

I converted the mergedpeaks assay to a ChromatinAssay class(integrated@assays$mergedpeaks) [1] "ChromatinAssay" attr(,"package") [1] "Signac"

however i am still getting an error when using coverageplot()

CoveragePlot(

• object = integrated,
• region = "chr2-87011729-87035519"
• ) Error in colnames<-(*tmp*, value = start(x = region):end(x = region)) : attempt to set 'colnames' on an object with less than two dimensions

ensuring defaultassay is mergedpeaks gives the same error

DefaultAssay(integrated) <- 'mergedpeaks' CoveragePlot(

• object = integrated,
• region = "chr2-87011729-87035519"
• ) Error in colnames<-(*tmp*, value = start(x = region):end(x = region)) : attempt to set 'colnames' on an object with less than two dimensions

[timoast]

This might be due to no Tn5 integration events are being found in the requested region. Can you confirm:

1. The coordinates are correct (eg, chromosome names in the fragment file start with "chr")
2. There are fragments present in the region requested (you can check this using tabix with the fragment file)
3. The fragment file information is set correctly for the "mergedpeaks" assay. You can extract the fragment information using Fragments(integrated[["mergedpeaks"]]), and check the fragment file path for each fragment object using the GetFragmentData function
4. The cell names in the fragment objects match the cell names in the fragment file. You can check the cell names in the fragment object using the GetFragmentData function with slot="cells".

[afcmalone]

Thank you. i will check these points and let you know if i run into more problems.

[achieng33]

I'm running into the same problem, and I've realized that my combined peaks assay doesn't contain a fragment file (it gives me a null list), so I went back to my script for when I created a chromatin assay object with the new file. Here's my object:

`counts_F1 <- Read10X_h5(filename = "/Users/Nils/OneDrive/MaryAnne R/P0 mice ATAC-seq/F1/filtered_peak_bc_matrix.h5") #peak/cell matrix; rows represent peaks (regions of open chromatin) and each value in the matrix is the number of Tn5 integration sites for each cell that map within each peak

metadata_F1 <- read.csv( file = "/Users/Nils/OneDrive/MaryAnne R/P0 mice ATAC-seq/F1/singlecell.csv", header = TRUE, row.names = 1 #contains the barcodes i.e. cell identities )

chrom_assay_F1 <- CreateChromatinAssay( counts = counts_F1, sep = c(":", "-"), genome = 'mm10', fragments = '/Users/Nils/OneDrive/MaryAnne R/P0 mice ATAC-seq/F1/fragments.tsv.gz', min.cells = 1 )

P0_F1 <- CreateSeuratObject( counts = chrom_assay_F1, assay = 'peaks', project = 'ATAC', meta.data = metadata_F1 )

P0_F1`

This is one of 4 objects I need to integrate. I used this code to form a common peak set:

`# read in peak sets - the BED file with the peaks was downloaded

peaks.F1 <- read.table( file = "/Users/Nils/OneDrive/MaryAnne R/P0 mice ATAC-seq/F1/peaks.bed", col.names = c("chr", "start", "end") #specifies names of variables ) peaks.F2 <- read.table( file = "/Users/Nils/OneDrive/MaryAnne R/P0 mice ATAC-seq/F2/peaks.bed", col.names = c("chr", "start", "end") ) peaks.M1 <- read.table( file = "/Users/Nils/OneDrive/MaryAnne R/P0 mice ATAC-seq/M1/peaks.bed", col.names = c("chr", "start", "end") ) peaks.M2 <- read.table( file = "/Users/Nils/OneDrive/MaryAnne R/P0 mice ATAC-seq/M2/peaks.bed", col.names = c("chr", "start", "end") )

convert to genomic ranges - this function takes a data-frame-like object as input and tries to automatically find the columns that describe genomic ranges. It returns them as a GRanges object

gr.F1 <- makeGRangesFromDataFrame(peaks.F1) gr.F2 <- makeGRangesFromDataFrame(peaks.F2) gr.M1 <- makeGRangesFromDataFrame(peaks.M1) gr.M2 <- makeGRangesFromDataFrame(peaks.M2)

Create a unified set of peaks to quantify in each dataset

combined.peaks <- reduce(x = c(gr.F1, gr.F2, gr.M1, gr.M2))

Filter out bad peaks based on length

peakwidths <- width(combined.peaks)

combined.peaks <- combined.peaks[peakwidths < 10000 & peakwidths > 20]

combined.peaks `

so then I quantify the merged peak set and added the assay:

`# count fragments per cell overlapping the common set of peaks in P0_F1

com_peaks_F1 <- FeatureMatrix( fragments = Fragments(P0_F1), features = combined.peaks, cells = colnames(P0_F1) )

create a new assay and add it to the P0_F1 dataset

P0_F1[['ComPeaks']] <- CreateChromatinAssay( counts = com_peaks_F1, min.cells = 1, ranges = combined.peaks, genome = 'mm10' )

DefaultAssay(P0_F1) <- 'ComPeaks' P0_F1 <- RunTFIDF(P0_F1) `

however, this is what it shows me when I try to extract the fragment information:

`> DefaultAssay(P0_F1) <- 'peaks'

Fragments(P0_F1) [[1]] A Fragment object for 7495 cells

DefaultAssay(P0_F1) <- 'ComPeaks' Fragments(P0_F1) list()`

would this be why it's giving me the error message? I'm not sure how I would go about fixing this. Message is - (Error in colnames<-(tmp, value = start(x = region):end(x = region)) : attempt to set 'colnames' on an object with less than two dimensions)

[kinali]

Hi, I am getting the same error when I am trying to use Coverageplot with my integrated data (RNA+ATAC). As you suggested I checked the fragments and I got the same error too:

Fragments(integrated[["peaks"]])

Error in UseMethod(generic = "Fragments", object = object) : no applicable method for 'Fragments' applied to an object of class "Assay"

How can I fix this?

[timoast]

@kinali it looks like you're trying to run it on a standard Assay rather than a ChromatinAssay. You probably need to change the default assay in the Seurat object

[kinali]

Hi @timoast thank you for your response! Actually, I am trying to run it on ChromatinAssay (peaks is the name of the ChromatinAssay in here). I created ChromatinAssay before the integration like this:

chrom_assay <- CreateChromatinAssay( counts = counts, sep = c(":", "-"), genome = 'hg38', fragments = './fragments.tsv.gz', min.cells = 10, min.features = 200 ) s1 <- CreateSeuratObject(counts = chrom_assay, assay = "peaks", meta.data = metadata, project = "s1") Before integration it is looking like Chromatin assay: s1@assays$peaks ChromatinAssay data with 123082 features for 9360 cells Variable features: 123082 Genome: hg38 Annotation present: TRUE Motifs present: FALSE Fragment files: 1

After integration it changes to Chromatinassay to assay.all is the name of the integrated object and I am getting that error:

`> all[["peaks"]] Assay data with 123082 features for 26135 cells First 10 features: chr1-181281-181662, chr1-629729-630172, chr1-633787-634275, chr1-778214-779318, chr1-817031-818204, chr1-818881-818960, chr1-826789-827920, chr1-867615-867878, chr1-869628-870243, chr1-904240-905757

DefaultAssay(all) <- "peaks" CoveragePlot(object = all,region = "chr10-129835283-129973053",features = "EBF3",peaks = TRUE ) Error in UseMethod(generic = "Fragments", object = object) : no applicable method for 'Fragments' applied to an object of class "Assay"`

My data is non-multiome, I tried to integrate separate 10X scRNA-seq and scATAC-seq objects. I have not used LinkPeaks. Could that be the reason?

[timoast]

@kinali that's a standard assay rather than a ChromatinAssay. If you need further support please open a new issue with your full code