在Matlab中探索基因表达数据

在Matlab中探索基因表达数据

本文利用Matlab及其生物信息学工具箱提供的函数识别差异表达基因并利用基因本体论确定差异表达基因的生物学功能。

引言

包含寡核苷酸或cDNA探针的微阵列可用来比较基因组尺度的基因表达谱，微阵列试验的重要目的在于确定不同条件下，如两种不同的肿瘤类型，是否存在统计显著的基因表达量的变化进而确定差异表达基因的生物学功能。

本文利用一个公共数据集来说明计算过程，这个数据集包括42个胚胎中枢神经系统肿瘤组织（CNS, Pomeroy et al. 2002），样本采用Affymetrix 公司出品的HuGeneFL基因芯片进行杂交。

这些CNS数据集（CEL文件）可在CNS实验网站获得，42个肿瘤样本包括10个10 个髓母细胞瘤, 10个横纹肌样脑膜瘤, 10个胶质瘤, 8个幕上原始神经外胚层肿瘤和4个正常人小脑，CNS原始数据集用鲁棒多芯片平均（RMA）和GC鲁棒多芯片平均（GCRMA）进行了预处理。可以采用t检验和假发现率（FDR）来检测不同肿瘤类型间差异表达的基因，还可以探索与显著上跳基因相关的基因本体论术语。

载入基因表达数据

用Load命令加载MAT文件cnsexpressiondata包含三个DataMatrix对象，expr_cns_rma, expr_cns_gcrma_mle, and expr_cns_gcrma_eb，分别储存用RMA和GCRMA（MLE和EB）预处理的基因表达值。

load cnsexpressiondata

在每个DataMatrix对象中，每行对应一个HuGeneFl芯片的探针集，每列对应于一个样本，行名是探针集的ID而列名为样本名，本文用expr_cns_gcrma_eb示例，当然也可以用其他对象。

调用get命令获取DataMatrix对象的特征。

get(expr_cns_gcrma_eb)

Name: 'CNS gene expression data'

RowNames: {7129x1 cell}

ColNames: {1x42 cell}

NRows: 7129

NCols: 42

NDims: 2

ElementClass: 'single'

确定DataMatrix对象expr_cns_gcrma_eb中的基因和样本的数目。

[nGenes, nSamples] = size(expr_cns_gcrma_eb)

nGenes =

7129

nSamples =

42

可以用基因符号来代替探针集的ID用于标记基因表达值，HuGeneFl芯片的基因符号在一个包含Java哈希表的MAT文件中。

load HuGeneFL_genesymbol_hashtable;

为hu6800genesymbol_hashtable变量创建一个基因表达值的基因符号的cell矩阵。

huGenes = cell(nGenes, 1);

for i =1:nGenes

huGenes{i} = hu6800genesymbol_hashtable.get(expr_cns_gcrma_eb.RowNames{i});

end

用DataMatrix的rownames方法将exprs_cns_gcrma_eb中的行名设成基因符号。

expr_cns_gcrma_eb = rownames(expr_cns_gcrma_eb, ':', huGenes);

基因表达数据的过滤

首先除去没有基因符号的表达数据，如标成'---'的空符号。

expr_cns_gcrma_eb('---', :) = [];

在这个研究中很多基因没有表达或在样本间变化很小，这些基因需要用非特异性过滤除去。用genelowvalfilter函数滤除绝对表达量值很低的基因。

[mask, expr_cns_gcrma_eb] = genelowvalfilter(expr_cns_gcrma_eb);

用genevarfilter函数滤除样本间方差很小的基因。

[mask, expr_cns_gcrma_eb] = genevarfilter(expr_cns_gcrma_eb);

确定过滤以后的基因数目。

nGenes = expr_cns_gcrma_eb.NRows

nGenes =

5758

识别差异基因表达

现在可以比较一下CNS髓母细胞瘤(MD)和非神经源恶性胶质瘤(Mglio)之间基因表达值的差异了。

从42个样本中提取10个MD和10个Mglio样本数据。

MDs = strncmp(expr_cns_gcrma_eb.ColNames,'Brain_MD', 8);

Mglios = strncmp(expr_cns_gcrma_eb.ColNames,'Brain_MGlio', 11);

MDData = expr_cns_gcrma_eb(:, MDs);

get(MDData)

Name: ''

RowNames: {5758x1 cell}

ColNames: {1x10 cell}

NRows: 5758

NCols: 10

NDims: 2

ElementClass: 'single'

MglioData = expr_cns_gcrma_eb(:, Mglios);

get(MglioData)

Name: ''

RowNames: {5758x1 cell}

ColNames: {1x10 cell}

NRows: 5758

NCols: 10

NDims: 2

ElementClass: 'single'