0

我对为什么这段代码不断抛出错误感到非常困惑。

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.

这非常令人困惑，有人可以帮助我理解吗？谢谢！

4

1 回答 1

0

问题是对象中的列的子集函数错误，例如：

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
于 2020-03-18T15:30:08.917 回答