文章目录

1、分类

数据：癌症基因图谱（ TCGA）数据库中 5 种不同类型癌症（包括BRCA, BLCA, LGG, LUAD 和 LUSC,每种为一个数据集）病例的 RNA 转录组（RNA-seq）数据。

关于数据集的下载在上一篇文章。下载完成后每一个病例都是一个.gz的压缩包。网上大部分关于这个数据库的处理都是用的R语言，自己想尝试一下用Python处理一下。

将所有样例压缩包放在同一文件下解压，将3329个TCGA样本通过.tsv文件进行整理，分成BRCA, BLCA, LGG, LUAD 和 LUSC五个数据集，其中，LGG有529个，LUAD有594个，LUSC有551个，BLCA有433个，BRCA有1222个。

import pandas as pdimport numpy as npimport osimport shutil#1.读取.tsv数据文件，将病例的文件名与癌肿名一一对应sample_s= pd.read_csv('gdc_sample_sheet.-03-16.tsv', sep='\t')file_dir = "D:\\python\\data\\data_yo\\" + sample_s.iloc[1, 1]#D:\\python\\data\\data_yo\\是.tsv文件保存的路径，sample_sum.iloc[1, 1]是文件保存的名字# for root, dirs, files in os.walk(file_dir):#for file in files:## file_obj=open(file_dir + "\\" + file)#读取文件内容#开始更改文件的路径for i in range(3329):if sample_s.iloc[i, 4] == 'TCGA-LGG':filename = sample_s.iloc[i, 1]file_dir = "D:\\python\\data\\data_yo\\" + filenameif os.path.exists(file_dir):#判断是否存在，存在就更改路径shutil.move(file_dir, "D:/python/data/LGG/")

2、数据整理

分类完成后每一类中都存在许多独立病例文件，而且文件名杂乱，不方便后续的数据分析，于是在处理分析数据前先对数据进行前处理。

2.1 换ID名

为了便于区分病例样本，将病例的ID和文件对应更改，这里需要用到.json文件。

这里借鉴了其他作者的思路，完整代码如下：

import osimport json#打开json文件夹；file = open('D:\\python\\data\\metadata.cart.-03-20.json', encoding='utf-8')#读取json文件夹；json_precess = json.loads(str(file.read()))#创建空字典；dict ={}for i in json_precess:print(i['file_name'])print(i['associated_entities'][0]['entity_submitter_id'])dict[str(i['file_name']).strip('.gz')] =i['associated_entities'][0]['entity_submitter_id']print(dict)#mainfest文件路径；path ='D:\\python\\data\\gdc_download_0316_024702.575784'filelist = os.listdir(path)#mainfest子路径下的所有文件列表；for file_one in filelist:file = path + '/' +file_oneprint(file)list = os.listdir(file)[0]print(list)if '.gz' in list:olddir = file +'/' +list#原来文件名newdir = file + '/' + dict[list.split('.gz')[0]] +'.gz'#新的文件名os.rename(olddir,newdir)#重新命名

2.2矩阵整理

遍历文件夹下的所有病例文件，因为每一个病例在相同位置对应的RNA基因是一样的，这样就方便读取并取出数据。将RNA基因名放在第一列。最后得到5个矩阵，每个代表一个类。

import gzipimport jsonimport osimport sysimport reimport shutilfrom collections import OrderedDictfrom datetime import timefrom multiprocessing import Poolimport globimport asyncioimport shutilfile_dir1=('D:\\python\\data\\data_BLCA')file_dir2=('D:\\python\\data\\data_BRCA')file_dir3=('D:\\python\\data\\data_LGG')file_dir4=('D:\\python\\data\\data_LUAD')file_dir5=('D:\\python\\data\\data_LUSC')#读取每个文件夹下解压后的.txt文件内容，将所有数据整合到一个矩阵中def readfile(file):# file_obj=open(file)file_obj=open(file_dir5 + "\\" + file)try:while True:line = file_obj.readline().strip("\n")if not line:breakarry = line.split("\t")if arry[0] in tcga_matrix:tcga_matrix[arry[0]].append(arry[1])else:tcga_matrix[arry[0]] = [arry[1]]finally:file_obj.close()if __name__ =="__main__":tcga_matrix=OrderedDict()for root, dirs, files in os.walk(file_dir5):for file in files:# file_obj=open(file_dir + "\\" + file)#读取文件内容readfile(file)with open('LUSC_matrix.txt','w') as f_w:for k,v in tcga_matrix.items():f_w.write(k+'\t'+'\t'.join(v[0:len(v)])+'\n')# f_w.write(k+'\t'+'\t'.join(v[0:len(v)])+'\n')f_w.close()

LGG部分数据：

（a）对于每一个数据集，在 6 万多个 RNA 中，如 果某 RNA 在一半以上的病例中没有表达（数值为 0），则将其丢弃。

import matplotlib.pyplot as pltimport numpy as npimport pandas as pdfrom numpy import aroundfrom pandas import Series,DataFrameimport palettableimport random#将超过一半数据为0的RNA删除掉，剩下3万左右数据file_dir1=("D:\\python\\data\\BLCA_matrix.txt")file_dir2=("D:\\python\\data\\BRCA_matrix.txt")file_dir3=("D:\\python\\data\\LGG_matrix.txt")file_dir4=("D:\\python\\data\\LUSC_matrix.txt")file_dir5=("D:\\python\\data\\LUAD_matrix.txt")aaa=0#计算有多少RNA没被删掉## #新文件的地址f1 = open('LUAD.txt', "w") # 以写入的形式打开txt文件with open(file_dir5, "r") as f:#读取.txt每一行for line in f.readlines():# print(line)#每个数据以\t隔开sub_str = line.split('\t')# print(sub_str[0])# print(sub_str[1])## print(type(sub_str))count = 0#每一个RNA基因含0的个数num_list_new = [float(x) for x in sub_str[1:]]#将list-str转为list-float，不包含第一个（RNA名称）for i in num_list_new[0:]:#每一行数据if i == 0:count += 1# print(count)# print(len(num_list_new))length = len(num_list_new)# print(around(5/10))if around(count / length) == 0:#少于一半为0则转移到新的文件aaa+=1print(aaa)# ','.join(v[0:len(v)]) + '\n'# print(','.join(sub_str[0:len(sub_str)])) # 将修改后的文本内容写入f1.writelines('\t'.join(sub_str[0:len(sub_str)])) # 将修改后的文本内容写入f1.close() # 关闭文件

（b）对 RNA 的表达值进行对数标准化：𝑎 → 𝑙𝑜𝑔10(𝑎 + 1)。

a=math.log10(a + 1)