Ошибка в eval (expr = expr): объект 'nFeature_RNA_lower' не найден

Question

Я очень озадачен тем, почему этот код продолжает выдавать ошибку.

qc_and_normalize <- function(seurat_object, coeff=3) {

  nCount_RNA_upper <- median(seurat_object@meta.data$nCount_RNA) + coeff*mad(seurat_object@meta.data$nCount_RNA)

  nFeature_RNA_upper <- median(seurat_object@meta.data$nFeature_RNA) +   coeff*mad(seurat_object@meta.data$nFeature_RNA)

  nCount_RNA_lower <- median(seurat_object@meta.data$nCount_RNA) - coeff*mad(seurat_object@meta.data$nCount_RNA)

  nFeature_RNA_lower <- median(seurat_object@meta.data$nFeature_RNA) - coeff*mad(seurat_object@meta.data$nFeature_RNA)


  dims <- seurat_object@assays$RNA@counts@Dim

  dims_lower_cut <- sum((seurat_object$nCount_RNA < nCount_RNA_lower) & (seurat_object$nFeature_RNA < nFeature_RNA_lower))
  dims_upper_cut <- sum((seurat_object$nCount_RNA > nCount_RNA_upper) & (seurat_object$nFeature_RNA > nFeature_RNA_upper))
  print(paste0('removing ', round((dims_lower_cut + dims_upper_cut) / dims[2] * 100, 2), '% of samples'))

  # QC
  # remove "bad cells"
  seurat_object <- subset(seurat_object, ( nFeature_RNA > nFeature_RNA_lower | nCount_RNA > nCount_RNA_lower) 
                          & (nCount_RNA < nCount_RNA_upper | nFeature_RNA < nFeature_RNA_upper))


  # normalize data
  # by default, normalization.method = "LogNormalize", scale.factor = 10000
  seurat_object <- NormalizeData(seurat_object)

  return(seurat_object)
}

Функция не так важна, просто когда она достигает линии

seurat_object <- subset(seurat_object, ( nFeature_RNA > nFeature_RNA_lower | nCount_RNA > nCount_RNA_lower) 
                          & (nCount_RNA < nCount_RNA_upper | nFeature_RNA < nFeature_RNA_upper))

выдает ошибку

Error in eval(expr = expr) : object 'nFeature_RNA_lower' not found

Почему она не может найти nFeature_RNA_lower когда это было определено в функции? Когда я пытаюсь напечатать его над строкой ошибки, он печатает его без проблем. Все работает также, когда я запускаю код построчно в глобальной среде, но когда я превращаю его в функцию, он внезапно не может найти nFeature_RNA_lower.

Это очень запутанно, может кто-нибудь помочь мне понять? Спасибо!

Answer 1

Проблема в подмножестве функций ошибок, которые для столбца в объекте, например:

neu = function(x){
 cutoff = 3
 subset(x,nFeature_RNA>cutoff)
 }

neu(pbmc_small)
Error in eval(expr = expr) : object 'cutoff' not found

Вы можете подмножество по матрице:

neu = function(x){
cutoff = 3
x[,x@meta.data$nFeature_RNA>cutoff]
}
neu(pbmc_small)
An object of class Seurat 
230 features across 80 samples within 1 assay 
Active assay: RNA (230 features)
 2 dimensional reductions calculated: pca, tsne

Глядя на функцию у вас есть, мы также можем немного упростить это, сначала мы определяем сэмплы, которые хотим исключить:

foo = function(features,coeff){
  m = apply(features,2,function(i){
    abs(i - median(i)) > coeff*mad(i)
  })
  rowMeans(m) == 1
}

Затем необходимо использовать это в основной функции:

qc_and_normalize <- function(seurat_object, coeff=2) {
  exclude = foo(seurat_object@meta.data[,c("nCount_RNA","nFeature_RNA")],coeff)
  seurat_object <- seurat_object[,!exclude]
  seurat_object <- NormalizeData(seurat_object)
  return(seurat_object)
}

Мы можем проверить это:

qc_and_normalize(pbmc_small,2)
Performing log-normalization
0%   10   20   30   40   50   60   70   80   90   100%
[----|----|----|----|----|----|----|----|----|----|
**************************************************|
An object of class Seurat 
230 features across 77 samples within 1 assay 
Active assay: RNA (230 features)
 2 dimensional reductions calculated: pca, tsne

Проверьте это, используя свой код:

seurat_object = pbmc_small
nCount_RNA_upper <- median(seurat_object@meta.data$nCount_RNA) + coeff*mad(seurat_object@meta.data$nCount_RNA)
nFeature_RNA_upper <- median(seurat_object@meta.data$nFeature_RNA) +   coeff*mad(seurat_object@meta.data$nFeature_RNA)
nCount_RNA_lower <- median(seurat_object@meta.data$nCount_RNA) - coeff*mad(seurat_object@meta.data$nCount_RNA)
nFeature_RNA_lower <- median(seurat_object@meta.data$nFeature_RNA) - coeff*mad(seurat_object@meta.data$nFeature_RNA)

  dims_lower_cut <- sum((seurat_object$nCount_RNA < nCount_RNA_lower) & (seurat_object$nFeature_RNA < nFeature_RNA_lower))
  dims_upper_cut <- sum((seurat_object$nCount_RNA > nCount_RNA_upper) & (seurat_object$nFeature_RNA > nFeature_RNA_upper))

dims_lower_cut
[1] 0
> dims_upper_cut
[1] 3