python - 在 python 中用大文件绘制集群热图？

Question

我有三个文件，每个文件大小约为 3Gb，包含如下基因组映射信息：

#RefName    Pos Coverage
BGC0000001_59320bp  0   0
BGC0000001_59320bp  1   0
BGC0000001_59320bp  2   0
BGC0000001_59320bp  3   0
BGC0000001_59320bp  4   0
BGC0000001_59320bp  5   0
BGC0000001_59320bp  6   0
BGC0000001_59320bp  7   0
BGC0000001_59320bp  8   0
BGC0000001_59320bp  9   0
BGC0000001_59320bp  10  0
BGC0000001_59320bp  11  0
BGC0000001_59320bp  12  0
BGC0000001_59320bp  13  0
BGC0000001_59320bp  14  0
BGC0000001_59320bp  15  0
BGC0000001_59320bp  16  0
BGC0000001_59320bp  17  0
BGC0000001_59320bp  18  0
BGC0000001_59320bp  19  0
BGC0000001_59320bp  20  0
BGC0000001_59320bp  21  0
BGC0000001_59320bp  22  0
BGC0000001_59320bp  23  0
BGC0000001_59320bp  24  0
BGC0000001_59320bp  25  0
BGC0000001_59320bp  26  0
BGC0000001_59320bp  27  0
BGC0000001_59320bp  28  0
BGC0000001_59320bp  29  0
BGC0000001_59320bp  30  0
BGC0000001_59320bp  31  0
BGC0000001_59320bp  32  0
BGC0000001_59320bp  33  0
BGC0000001_59320bp  34  0
BGC0000001_59320bp  35  0

.
.
.
.

第一列 = 基因名称；第 2 列 = 基因每个碱基的位置或基因座；第三列 = 每个基地或位置的覆盖范围。

例如：第一个基因'BGC0000001_59320bp'是59320bp长，所以这个基因会有59320行，然后下一个基因出现。

我想制作一个像这样的集群热图：

y 轴是每个基因，例如 BGC0000001,BGC0000002,.... x 轴是 binned 基因位置。因为每个基因都有不同的长度，我会将它们放入相同数量的箱中，比如 100 个箱。这使得 x 轴变为 bin1, bin2, bin3, ..., bin100。例如，BGC0000001 = 59320 bp，BGC0000002= 10000bp，我将每个基因分成 100 个 bin，所以对于 BGC0000001，bin1 = 593bp，bin2 = 593bp，...bin100= 595bp；对于 BGC0000002，bin1 = 100bp，bin2 = 100bp ...bin100 = 100bp。我会将覆盖范围的总和作为绘图值放入每个 bin 中。例如：对于BGC0000001，bin1 = sum(前593bp的覆盖率)，bin2 = sum(第二个593bp的覆盖率) ...

最后，我将使用分箱数据进行聚类热图绘制。我有三个示例文件来制作三个图。因为数据量很大，我正在尝试下面的代码，确实不干净，而且需要很长时间：

import pandas as pd
import matplotlib.pyplot as plt
import sys
import os
from scipy.stats import binned_statistic

import seaborn as sns



def getFiles(filePath):
    coverageFile = []    

    for file in os.listdir(filePath):
        if file.endswith(".coverage"):
            coverageFile.append(file)
        else:
            pass          

    print(coverageFile)

    return coverageFile



def binCov(df,name,nBins):
    genomeMap = df[df.node == name].coverage
    bin_sums, bin_edges, binnumber = binned_statistic (range(len(genomeMap)), genomeMap, 
                                                        statistic = 'sum', bins = nBins)
    bin_width = (bin_edges[1] - bin_edges[0])
    bin_centers = bin_edges[1:] - bin_width*0.5

    df1 = pd.DataFrame({'covSum':bin_sums})

    df1.insert(1,'gene',name)
    df1.insert(2,'loci',range(1, len(df1)+1))
    df1.insert(2,'covNorm',df1.covSum)
    df1.loc[df1.covSum >5000, 'covNorm'] = 5000
    df1.drop('covSum', axis =1, inplace = True)

    return df1


files = getFiles(os.getcwd())

df_nac = pd.read_csv(files[0], sep = '\t')
df_nac.columns = ['node', 'loci', 'coverage']


df_ks = pd.read_csv(files[1], sep = '\t')
df_ks.columns = ['node', 'loci', 'coverage']

df_meta = pd.read_csv(files[2], sep = '\t')
df_meta.columns = ['node', 'loci', 'coverage']

cluster_names = df_nac.node.tolist()

cluster_nac = {}
cluster_ks = {}
cluster_meta = {}


for name in cluster_names:
    cluster_nac[name] = binCov(df_nac, name, 100)
    cluster_ks[name] = binCov(df_ks, name, 100)
    cluster_meta[name] = binCov(df_meta, name, 100)

cluster_nac_all = pd.concat(cluster_nac.values(), ignore_index=True)
cluster_ks_all = pd.concat(cluster_ks.values(), ignore_index=True)
cluster_meta_all = pd.concat(cluster_meta.values(), ignore_index=True)


nac_p = cluster_nac_all.pivot("gene", "loci", "covNorm")
ks_p = cluster_ks_all.pivot("gene", "loci", "covNorm")
meta_p = cluster_meta_all.pivot("gene", "loci", "covNorm")



sns.set_style("darkgrid")

plt.figure(figsize = (5,10))
# cmap="Greens", , linecolor='yellow', alpha = 0.9,cmap="YlGnBu"
g1 = sns.heatmap(meta_p, xticklabels=True, yticklabels=True, cmap="BuPu", cbar = False)

g1.set_xlabel('Bases (binned)')
g1.set_ylabel('')

g1.set_xticks([])
g1.set_yticks([])

plt.savefig('meta_sumCov.png',dpi = 300, bbox_inches = 'tight')

plt.figure(figsize = (5,10))
# cmap="Greens", , linecolor='yellow', alpha = 0.9,cmap="YlGnBu"
g2 = sns.heatmap(nac_p, xticklabels=True, yticklabels=True, cmap="BuPu", cbar = False)

g2.set_xlabel('Bases (binned)')
g2.set_ylabel('')

g2.set_xticks([])
g2.set_yticks([])

plt.savefig('nac_sumCov.png',dpi = 300, bbox_inches = 'tight')

plt.figure(figsize = (5,10))
# cmap="Greens", , linecolor='yellow', alpha = 0.9,cmap="YlGnBu"
g3 = sns.heatmap(ks_p, xticklabels=True, yticklabels=True, cmap="BuPu", 
    cbar_kws = dict(use_gridspec=False,location="right"))

g3.set_xlabel('Bases (binned)')
g3.set_ylabel('')

g3.set_xticks([])
g3.set_yticks([])

plt.savefig('ks_sumCov.png',dpi = 300, bbox_inches = 'tight')

我让它运行了 10 个小时，但它仍然运行，尽管它没有显示任何错误。我不认为 Pandas 对大文件有好处，但我不熟悉其他更快的工具。请帮我。

谢谢你。

标签： pythonheatmap