scanpy赋值问题

2023-11-03 01:20
文章标签 问题 赋值 scanpy

本文主要是介绍scanpy赋值问题,希望对大家解决编程问题提供一定的参考价值,需要的开发者们随着小编来一起学习吧!

今天发现一个很奇怪的bug


import numpy as np
import pandas as pd
import anndata as ad
from scipy.sparse import csr_matrix
print(ad.__version__)counts = csr_matrix(np.random.poisson(1, size=(100, 2000)), dtype=np.float32)
adata1 = ad.AnnData(counts)
print(adata1)def f(adata):adata = adata[:,0:1] # print(adata.shape)f(adata1)
print(adata1.shape)

结果如下
在这里插入图片描述
可以看到在函数中,这个adata的结果是变化了,但是并没有改变外部adata的值


import numpy as np
import pandas as pd
import anndata as ad
from scipy.sparse import csr_matrix
print(ad.__version__)counts = csr_matrix(np.random.poisson(1, size=(100, 2000)), dtype=np.float32)
adata1 = ad.AnnData(counts)
print(adata1.X[0:2,0:10])def f(adata):adata = adata[:,0:1] # print(adata.shape)f(adata1)
print(adata1.shape)
print(adata1.X[0:2,0:10])

在这里插入图片描述

但是如果一开始我不在函数中操作,而是主程序中,这个结果


import numpy as np
import pandas as pd
import anndata as ad
from scipy.sparse import csr_matrix
print(ad.__version__)counts = csr_matrix(np.random.poisson(1, size=(100, 2000)), dtype=np.float32)
adata1 = ad.AnnData(counts)
print(adata1.X.shape)adata1 = adata1[:,0:1]
print(adata1.shape)

结果如下
在这里插入图片描述
这个现象只能解释为adata= adata1[:,0:1]是一个复制的行为,只不过同名了,所以adata的饮用变了,如果

adata2 = adata1[:,0:1],

可以想象,这个结果不会对adata1结果有影响

这仅仅是一个简简单单的例子,下面有一个更奇怪的测试

import scanpy as sc adata= sc.read("/Users/yxk/Desktop/test_dataset/pbmc/pbmc.h5ad")
adata.obs["BATCH"] = adata.obs["batch"].copy()
adata.obs["label"]=adata.obs["celltype"].astype("category").cat.codes 
n_classes= len(adata.obs["label"].value_counts())
print(adata)adata1= adata[adata.obs["batch"]=="pbmc_3p"].copy()
adata2= adata[adata.obs["batch"]=="pbmc_5p"].copy()
#print(adata1.X)
#print(adata2.X)## 如果用这种方式,我的结果是这样的
def preprocessNew(adata_A_input, ):'''Performing preprocess for a pair of datasets.To integrate multiple datasets, use function preprocess_multiple_anndata in utils.py'''adata_A = adata_A_inputprint("Finding highly variable genes...")#sc.pp.highly_variable_genes(adata_A, flavor='seurat_v3', n_top_genes=2000)#hvg_A = adata_A.var[adata_A.var.highly_variable == True].sort_values(by="highly_variable_rank").indexprint("Normalizing and scaling...")sc.pp.normalize_total(adata_A, target_sum=1e4)sc.pp.log1p(adata_A)sc.pp.highly_variable_genes(adata_A,n_top_genes=2000)hvg_A = list(adata1.var_names[adata1.var.highly_variable])adata_A = adata_A[:, hvg_A]sc.pp.scale(adata_A, max_value=10)print(adata_A.X[0:1,0:100])print(adata_A.X.shape)# 为啥这些结果是这样的preprocessNew(adata1)
print(adata1.X.shape)

在这里插入图片描述可以看到adata的结果是没有改变的,还是33694维,但是我在函数中,明明是选择了高变基因的

但是如果采用下面的代码

import scanpy as sc adata= sc.read("/Users/yxk/Desktop/test_dataset/pbmc/pbmc.h5ad")
adata.obs["BATCH"] = adata.obs["batch"].copy()
adata.obs["label"]=adata.obs["celltype"].astype("category").cat.codes 
n_classes= len(adata.obs["label"].value_counts())
print(adata)adata1= adata[adata.obs["batch"]=="pbmc_3p"].copy()
adata2= adata[adata.obs["batch"]=="pbmc_5p"].copy()
#print(adata1.X)
#print(adata2.X)def preprocessNew(adata_A_input, ):'''Performing preprocess for a pair of datasets.To integrate multiple datasets, use function preprocess_multiple_anndata in utils.py'''adata_A = adata_A_inputprint("Finding highly variable genes...")#sc.pp.highly_variable_genes(adata_A, flavor='seurat_v3', n_top_genes=2000)#hvg_A = adata_A.var[adata_A.var.highly_variable == True].sort_values(by="highly_variable_rank").indexprint("Normalizing and scaling...")sc.pp.normalize_total(adata_A, target_sum=1e4)sc.pp.log1p(adata_A)sc.pp.highly_variable_genes(adata_A,n_top_genes=2000,subset=True)#adata_A = adata_A[:, hvg_A]sc.pp.scale(adata_A, max_value=10)print(adata_A.X[0:1,0:100])
preprocessNew(adata1)
print(adata1.X.shape)
print(adata1.X[0:1,0:100])
## 但是线则这个问题为啥不是

结果如下
在这里插入图片描述
这里可以看到,我最终的adata1的维度是改变了,这里需要注意

这里使用
sc.pp.highly_variable_genes(adata1,n_top_genes=2000,subset=True),就是对adata的引用改动了,最终导致最开始的atata出现了变化,反正最好还是用scanpy的内置函数了,一旦在函数里赋值就要注意局部对象的问题

import scanpy as sc adata= sc.read("/Users/yxk/Desktop/test_dataset/pbmc/pbmc.h5ad")
adata.obs["BATCH"] = adata.obs["batch"].copy()
adata.obs["label"]=adata.obs["celltype"].astype("category").cat.codes 
n_classes= len(adata.obs["label"].value_counts())
print(adata)
adata1= adata[adata.obs["batch"]=="pbmc_3p"].copy()
#adata2= adata[adata.obs["batch"]=="pbmc_5p"].copy()
#print(adata1.X)
#print(adata2.X)
print("Normalizing and scaling...")
sc.pp.normalize_total(adata1, target_sum=1e4)
sc.pp.log1p(adata1)
sc.pp.highly_variable_genes(adata1,n_top_genes=2000,subset=True)
sc.pp.scale(adata1, max_value=10)
print(adata1.X[0:1,0:100])
print(adata1.X.shape)
print(adata1.X[0:1,0:100])
## 但是线则这个问题为啥不是

如果采用了preprocessNew的函数,那么本质上只对adata做了如下变化

import scanpy as sc adata= sc.read("/Users/yxk/Desktop/test_dataset/pbmc/pbmc.h5ad")
adata.obs["BATCH"] = adata.obs["batch"].copy()
adata.obs["label"]=adata.obs["celltype"].astype("category").cat.codes 
n_classes= len(adata.obs["label"].value_counts())
print(adata)adata1= adata[adata.obs["batch"]=="pbmc_3p"].copy()
#adata2= adata[adata.obs["batch"]=="pbmc_5p"].copy()
#print(adata1.X)
#print(adata2.X)## 如果用这种方式,我的结果是这样的
def preprocessNew(adata_A_input, ):'''Performing preprocess for a pair of datasets.To integrate multiple datasets, use function preprocess_multiple_anndata in utils.py'''adata_A = adata_A_inputprint("Finding highly variable genes...")#sc.pp.highly_variable_genes(adata_A, flavor='seurat_v3', n_top_genes=2000)#hvg_A = adata_A.var[adata_A.var.highly_variable == True].sort_values(by="highly_variable_rank").indexprint("Normalizing and scaling...")sc.pp.normalize_total(adata_A, target_sum=1e4)sc.pp.log1p(adata_A)sc.pp.highly_variable_genes(adata_A,n_top_genes=2000)hvg_A = list(adata1.var_names[adata1.var.highly_variable])adata_A = adata_A[:, hvg_A]sc.pp.scale(adata_A, max_value=10)print(adata_A.X[0:1,0:100])print(adata_A.X.shape)# 为啥这些结果是这样的preprocessNew(adata1)
print(adata1.X.shape)
print(adata1.X[0:1,0:100])

结果如下
在这里插入图片描述reproduce result

import scanpy as sc adata= sc.read("/Users/yxk/Desktop/test_dataset/pbmc/pbmc.h5ad")
adata.obs["BATCH"] = adata.obs["batch"].copy()
adata.obs["label"]=adata.obs["celltype"].astype("category").cat.codes 
n_classes= len(adata.obs["label"].value_counts())
print(adata)adata2= adata[adata.obs["batch"]=="pbmc_3p"].copy()#print(adata1.X)
#print(adata2.X)## 如果用这种方式,我的结果是这样的print("Normalizing and scaling...")
sc.pp.normalize_total(adata2, target_sum=1e4)
sc.pp.log1p(adata2) # 真正对adata1只有这么多的操作# 为啥这些结果是这样的
print(adata2.X.shape)
print(adata2.X[0:1,0:100])

在这里插入图片描述

from sklearn.metrics import mean_squared_error
mean_squared_error(adata1.X.toarray(),adata2.X.toarray())

结果如下
在这里插入图片描述

这篇关于scanpy赋值问题的文章就介绍到这儿,希望我们推荐的文章对编程师们有所帮助!



http://www.chinasem.cn/article/334586

相关文章

MySQL 设置AUTO_INCREMENT 无效的问题解决

《MySQL设置AUTO_INCREMENT无效的问题解决》本文主要介绍了MySQL设置AUTO_INCREMENT无效的问题解决,文中通过示例代码介绍的非常详细,对大家的学习或者工作具有一定的参... 目录快速设置mysql的auto_increment参数一、修改 AUTO_INCREMENT 的值。

关于跨域无效的问题及解决(java后端方案)

《关于跨域无效的问题及解决(java后端方案)》:本文主要介绍关于跨域无效的问题及解决(java后端方案),具有很好的参考价值,希望对大家有所帮助,如有错误或未考虑完全的地方,望不吝赐教... 目录通用后端跨域方法1、@CrossOrigin 注解2、springboot2.0 实现WebMvcConfig

Go语言中泄漏缓冲区的问题解决

《Go语言中泄漏缓冲区的问题解决》缓冲区是一种常见的数据结构,常被用于在不同的并发单元之间传递数据,然而,若缓冲区使用不当,就可能引发泄漏缓冲区问题,本文就来介绍一下问题的解决,感兴趣的可以了解一下... 目录引言泄漏缓冲区的基本概念代码示例:泄漏缓冲区的产生项目场景:Web 服务器中的请求缓冲场景描述代码

Java死锁问题解决方案及示例详解

《Java死锁问题解决方案及示例详解》死锁是指两个或多个线程因争夺资源而相互等待,导致所有线程都无法继续执行的一种状态,本文给大家详细介绍了Java死锁问题解决方案详解及实践样例,需要的朋友可以参考下... 目录1、简述死锁的四个必要条件:2、死锁示例代码3、如何检测死锁?3.1 使用 jstack3.2

解决JSONField、JsonProperty不生效的问题

《解决JSONField、JsonProperty不生效的问题》:本文主要介绍解决JSONField、JsonProperty不生效的问题,具有很好的参考价值,希望对大家有所帮助,如有错误或未考虑... 目录jsONField、JsonProperty不生效javascript问题排查总结JSONField

github打不开的问题分析及解决

《github打不开的问题分析及解决》:本文主要介绍github打不开的问题分析及解决,具有很好的参考价值,希望对大家有所帮助,如有错误或未考虑完全的地方,望不吝赐教... 目录一、找到github.com域名解析的ip地址二、找到github.global.ssl.fastly.net网址解析的ip地址三

MySQL版本问题导致项目无法启动问题的解决方案

《MySQL版本问题导致项目无法启动问题的解决方案》本文记录了一次因MySQL版本不一致导致项目启动失败的经历,详细解析了连接错误的原因,并提供了两种解决方案:调整连接字符串禁用SSL或统一MySQL... 目录本地项目启动报错报错原因:解决方案第一个:第二种:容器启动mysql的坑两种修改时区的方法:本地

springboot加载不到nacos配置中心的配置问题处理

《springboot加载不到nacos配置中心的配置问题处理》:本文主要介绍springboot加载不到nacos配置中心的配置问题处理,具有很好的参考价值,希望对大家有所帮助,如有错误或未考虑... 目录springboot加载不到nacos配置中心的配置两种可能Spring Boot 版本Nacos

Java中JSON格式反序列化为Map且保证存取顺序一致的问题

《Java中JSON格式反序列化为Map且保证存取顺序一致的问题》:本文主要介绍Java中JSON格式反序列化为Map且保证存取顺序一致的问题,具有很好的参考价值,希望对大家有所帮助,如有错误或未... 目录背景问题解决方法总结背景做项目涉及两个微服务之间传数据时,需要提供方将Map类型的数据序列化为co

如何解决Druid线程池Cause:java.sql.SQLRecoverableException:IO错误:Socket read timed out的问题

《如何解决Druid线程池Cause:java.sql.SQLRecoverableException:IO错误:Socketreadtimedout的问题》:本文主要介绍解决Druid线程... 目录异常信息触发场景找到版本发布更新的说明从版本更新信息可以看到该默认逻辑已经去除总结异常信息触发场景复