Python生物資訊學⑥繪制熱圖及火山圖

Python生物資訊學⑥繪制熱圖及火山圖

Python生物資訊學③提取差異基因

通過上Python生物資訊學③提取差異基因得到了該資料集的差異分析的兩個關鍵參數，1.差異倍數（foldchange）以及2.差異的P值。本篇目的是得到滿足差異倍數和差異P值的基因，同時進行可視化（包括差異分析常見的火山圖和熱圖）。

繪制火山圖

（1）第一步制作差異分析結果資料框

genearray = np.asarray(pvalue)
 
result = pd.DataFrame({'pvalue':genearray,'FoldChange':fold})
 
result['log(pvalue)'] = -np.log10(result['pvalue'])      

（2）第二步制作火山圖的準備工作

標明的差異基因标準是 I.差異倍數的絕對值大于1，II. 差異分析的P值小于0.05

result['sig'] = 'normal'
 
result['size']  =np.abs(result['FoldChange'])/10
 
result.loc[(result.FoldChange> 1 )&(result.pvalue < 0.05),'sig'] = 'up'
result.loc[(result.FoldChange< -1 )&(result.pvalue < 0.05),'sig'] = 'down'      

ax = sns.scatterplot(x="FoldChange", y="log(pvalue)",
                      hue='sig',
                      hue_order = ('down','normal','up'),
                      palette=("#377EB8","grey","#E41A1C"),
                      data=result)
ax.set_ylabel('-log(pvalue)',fontweight='bold')
ax.set_xlabel('FoldChange',fontweight='bold')      

篩選差異基因

fold_cutoff = 1
pvalue_cutoff = 0.05
 
filtered_ids = []
for i in range(0, number_of_genes):
    if (abs(fold[i]) >= fold_cutoff) and (pvalue[i] <= pvalue_cutoff):
        filtered_ids.append(i)
        
filtered = data2.iloc[filtered_ids,:]
print("Number of DE genes: ")
print(len(filtered.index))      

繪制熱圖

熱圖（heatmap）是生物學文章裡（尤其是RNA-seq相關論文）經常出現的圖檔。熱圖的用途一般有兩個。以RNA-seq為例，熱圖可以：1）直覺呈現多樣本多個基因的全局表達量變化；2）呈現多樣本或多基因表達量的聚類關系。熱圖一般使用顔色（例如紅綠的深淺）來展示多個樣本多個基因的表達量高低，既直覺又美觀。同時可以對樣本聚類或者對基因聚類。

sns.clustermap(filtered, cmap='RdYlGn_r', standard_scale = 0)      

（1）每一行為一個基因，每一列為一個sample。 （2）綠色代表相對低表達，紅色代表相對高表達。 （3）相對接近的樣本或者基因會聚類在一起，比如探針名為101695_at的基因在GSM130370相對高表達，而在GSM130366低表達。