Replies: 1 comment
-
Moved to #841 |
Beta Was this translation helpful? Give feedback.
0 replies
Sign up for free
to join this conversation on GitHub.
Already have an account?
Sign in to comment
-
Hello, I created an ATAC object and used the matrix cell and peaks to generate an input_cds for cicero followed by cicero_cds. Please see the error below
ERROR
Overlap QC metrics:
Cells per bin: 50
Maximum shared cells bin-bin: 44
Mean shared cells bin-bin: 0.933890030721768
Median shared cells bin-bin: 0
Error in stop_if_wrong_length("'seqnames'", ans_len) :
'seqnames' must have the length of the object to construct (1) or
length 1
Calls: ConnectionsToLinks ... makeGRangesFromDataFrame -> GRanges -> new_GRanges -> stop_if_wrong_length
Execution halted
CODE
library(Signac)
library(Seurat)
library(EnsDb.Hsapiens.v86)
library(BSgenome.Hsapiens.UCSC.hg38)
library(SeuratWrappers)
library(monocle3)
library(cicero)
library(SeuratWrappers)
library(ggplot2)
library(patchwork)
library(utils)
atac <- readRDS(paste(Dir_save,"atac84_label",sep=""))
DefaultAssay(atac) <- "peaks"
genome <- seqlengths(BSgenome.Hsapiens.UCSC.hg38)
genome <- genome[1:24]
genome.df <- data.frame("chr" = names(genome), "length" = genome)
#matrix
indata <- atac@assays$peaks@counts
indata@x[indata@x > 0] <- 1
#colnames
cellinfo <- as.data.frame(colnames(indata))
row.names(cellinfo) <- cellinfo$cellinfo
names(cellinfo) <- "cells"
#rownames
rows <- noquote(rownames(indata))
peakinfo <- do.call('rbind', strsplit(as.character(rows),'-',fixed=TRUE))
peakinfo <- data.frame(peakinfo)
colnames(peakinfo) <- c("chr", "bp1", "bp2")
peakinfo$site_name <- paste(peakinfo$chr, peakinfo$bp1, peakinfo$bp2, sep="_")
row.names(peakinfo) <- peakinfo$site_name
row.names(indata) <- row.names(peakinfo)
colnames(indata) <- row.names(cellinfo)
#input_cds
input_cds <- suppressWarnings(new_cell_data_set(indata,cell_metadata = cellinfo,gene_metadata = peakinfo))
input_cds <- monocle3::detect_genes(input_cds)
input_cds <- input_cds[Matrix::rowSums(exprs(input_cds)) != 0,]
print("cds object created")
#cicero_cds
input_cds <- estimate_size_factors(input_cds)
input_cds <- preprocess_cds(input_cds, method = "LSI")
input_cds <- reduce_dimension(input_cds, reduction_method = 'UMAP',
preprocess_method = "LSI")
umap_coords <- reducedDims(input_cds)$UMAP
cicero_cds <- make_cicero_cds(input_cds, reduced_coordinates = umap_coords)
print("cicero cds done")
#conn
conns <- run_cicero(cicero_cds,genome.df)
head(conns)
CCAN_assigns <- generate_ccans(conns)
head(CCAN_assigns)
links <- ConnectionsToLinks(conns = conns, ccans = CCAN_assigns)
head(links)
Links(atac) <- links
The last step with adding the Links gives me the error above. I would be thankful for your suggestions to resolve this.
Thank you!
Beta Was this translation helpful? Give feedback.
All reactions