1
$\begingroup$

I have one function that is for reading .h5ad files

ReadH5AD <- function(file) {
  # Load libraries
  library(anndata)
  library(reticulate)
  
  # Check if the file has the .h5ad extension
  if (tolower(tools::file_ext(file)) != "h5ad") {
    stop("Only h5ad file paths are accepted!")
  }
  
  # Import and read the h5ad file
  ad <- import("anndata", convert = FALSE)
  h5ad <- ad$read_h5ad(file)
  
  # Extract counts matrix as a dense matrix
  counts <- as.matrix(py_to_r(h5ad$X))
  
  # Extract cell names
  cell_names <- py_to_r(h5ad$obs$index)
  
  # Extract gene names
  gene_names <- py_to_r(h5ad$var$index)
  
  # Transpose the matrix
  counts <- t(counts)
  
  # Set cell names and gene names as row and column names
  colnames(counts) <- cell_names
  rownames(counts) <- gene_names
  
  return(counts)
}

Now there is another function which is used for file which are .mtx type along with barcode.tsv and gene.tsv.

In order to detect organism type this function is used

detectOrganism <- function(outdir,data){
  gene.loc <- file.path(outdir, list.files(path=outdir, pattern='*genes\\.tsv', recursive=FALSE)[1])
  features.loc <- file.path(outdir, list.files(path=outdir, pattern='*features\\.tsv\\.gz', recursive=FALSE)[1])
  tryCatch({
    # Flag to indicate if data is from CellRanger >= 3.0
    checkGeneFile <- file.exists(gene.loc)
    if (!checkGeneFile) {
      genes <- read.table(file = features.loc, header = FALSE, sep = '\t',nrow=1, row.names = NULL)
    }
    else{
      genes <- read.table(file = gene.loc, header = FALSE, sep = '\t',nrow=1, row.names=NULL)
    }        
    if(grepl("^ENSMUSG",genes[1,1])){data$organism="Mouse"}
    else if(grepl("^ENSG",genes[1,1])){data$organism="Human"}
    else{data$organism="NA";stop()}
  }, 
  interrupt = function(x){errorOccured<<-TRUE; speciesErrorNotification()},
  error = function(x){errorOccured<<-TRUE; speciesErrorNotification()}
  )
}

This works when the input type is matrix.mtx,gene.tsv along with barcode.tsv. As this can detect using grepl from the gene.tsv file.

I would like to implement the same for .h5ad input.

  gene_names <- py_to_r(h5ad$var$index)

This part can extract the gene_names from a given .h5ad input file but how do I pass this to the detectOrganism function.

I have here basically two file types .h5ad and .mtx.

So my simple goal is the add the same functionality when the input is .h5ad instead of .mtx file type.

Any suggestion or help would be really appreciated

$\endgroup$
2
  • 2
    $\begingroup$ Please edit your question and add some example input files that we can use to reproduce the problem. Also show us the expected output from those examples so we can know if we get it right. $\endgroup$
    – terdon
    Aug 12, 2023 at 12:08
  • $\begingroup$ do i need to add the h5ad file? $\endgroup$
    – kcm
    Aug 13, 2023 at 15:11

0

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Browse other questions tagged or ask your own question.