识别蛋白质的重要特征和分类
本例展示了如何对质谱数据进行分类，并使用一些统计工具来寻找潜在的疾病标志物和蛋白质组学模式诊断。
简介
血清蛋白质组学模式诊断可用于区分患有和没有疾病的患者的样本。使用表面增强激光解吸和电离(SELDI)蛋白质谱生成剖面模式。这项技术有可能改善癌症病理的临床诊断测试。目标是选择一组减少的测量或“特征”，可用于区分癌症患者和对照组患者。这些特征将是特定质量/电荷值下的离子强度水平。
数据进行预处理
本例中的卵巢癌数据集来自FDA-NCI临床蛋白质组学项目数据库．数据集使用WCX2蛋白阵列生成。数据集包括95个对照组和121个卵巢癌。该数据集的详细描述请参见[1]和[4]。
本例假设您已经拥有预处理过的数据OvarianCancerQAQCdataset.mat．但是，如果没有数据文件，可以按照示例中的步骤重新创建利用序贯并行计算的光谱批处理．
或者，您也可以运行脚本msseqprocessing.m．
目录(fullfile (matlabroot,“例子”，“bioinfo”，“主要”）)确保支持文件在搜索路径上。类型msseqprocessing
使用MSSEQPROCESSING脚本创建OvarianCancerQAQCdataset。mat(用于% CANCERDETECTDEMO)。在运行此文件之前，将变量% "repository"初始化为您放置质谱%文件的完整路径。例如:% % repository = 'F:/MassSpecRepository/OvarianCD_PostQAQC/';% % or % % repository = '/home/username/MassSpecRepository/OvarianCD_PostQAQC/';大约执行时间为18分钟(奔腾4,4ghz)。如果您有并行计算工具箱，请参阅BIODISTCOMPDEMO以了解如何加快分析速度。The MathWorks, Inc. repositoryC = [repository 'Cancer/'];repositoryN = [repository 'Normal/'];filesCancer = dir([repositoryC '*.txt']);NumberCancerDatasets = number (filesCancer); fprintf('Found %d Cancer mass-spectrograms.\n',NumberCancerDatasets) filesNormal = dir([repositoryN '*.txt']); NumberNormalDatasets = numel(filesNormal); fprintf('Found %d Control mass-spectrograms.\n',NumberNormalDatasets) files = [ strcat('Cancer/',{filesCancer.name}) ... strcat('Normal/',{filesNormal.name})]; N = numel(files); % total number of files fprintf('Total %d mass-spectrograms to process...\n',N) [MZ,Y] = msbatchprocessing(repository,files); disp('Finished; normalizing and saving to OvarianCancerQAQCdataset.mat.') Y = msnorm(MZ,Y,'QUANTILE',0.5,'LIMITS',[3500 11000],'MAX',50); grp = [repmat({'Cancer'},size(filesCancer));... repmat({'Normal'},size(filesNormal))]; save OvarianCancerQAQCdataset.mat Y MZ grp
上面列出的脚本和示例中的预处理步骤旨在说明一组有代表性的可能的预处理过程。使用不同的步骤或参数可能会导致本例的不同结果，并且可能会得到改进。
加载数据
一旦有了预处理数据，就可以将其加载到MATLAB中。
负载OvarianCancerQAQCdataset谁
名称大小字节类属性MZ 15000x1 120000双Y 15000x216 25920000双grp 216x1 25056 cell
有三个变量:MZ，Y，grp．MZ是质量/电荷矢量，Y是否所有216名患者(对照组和癌症患者)的强度值grp保存关于这些样本中哪些代表癌症患者，哪些代表正常患者的索引信息。
初始化一些将在整个示例中使用的变量。
N =数字(grp);样本数量%Cidx = strcmp(“癌症”、grp);癌症样本的逻辑指数载体Nidx = strcmp(“正常”、grp);正常样本的逻辑索引向量Cvec = find(Cidx);癌症样本的指数载体Nvec = find(Nidx);%正常样本的指数向量xAxisLabel =“质量/电荷(M / Z)”；% x标签用于绘图yAxisLabel =离子强度的；% y标签用于绘图
可视化一些样本
您可以将一些数据集绘制到图形窗口中，以直观地比较来自两个组的概要文件;在这个例子中，显示了癌症患者(蓝色)和对照组患者(绿色)的5个频谱图。
图;持有在；hC = plot(MZ,Y(:，Cvec(1:5))，“b”）;hN = plot(MZ,Y(:，Nvec(1:5))，‘g’）;包含(xAxisLabel);ylabel (yAxisLabel);轴([2000 12000 -5 60])图例([hN(1)，hC(1)]，{“控制”，卵巢癌的})标题(“多样本谱图”）
放大8500到8700 M/Z的区域，可以看到一些可能对数据分类有用的峰值。
轴([8450、8700、1、7])
另一种可视化整个数据集的方法是查看对照组和癌症样本的组平均信号。你可以画出组的平均值和每个组的信封。
mean_N = mean(Y(:，Nidx)，2);%组均值为对照样本max_N = max(Y(:，Nidx)，[]，2);%对照样品的顶部信封min_N = min(Y(:，Nidx)，[]，2);%对照样品的底部信封mean_C = mean(Y(:，Cidx)，2);癌症样本的%组平均值max_C = max(Y(:，Cidx)，[]，2);%对照样品的顶部信封min_C = min(Y(:，Cidx)，[]，2);%对照样品的底部信封图;持有在；hC = plot(MZ,mean_C，“b”）;hN = plot(MZ,mean_N，‘g’）;gC = plot(MZ，[max_C min_C]，“b——”）;gN = plot(MZ，[max_N min_N]，“g——”）;包含(xAxisLabel);ylabel (yAxisLabel);轴([8450、8700、1、7])传说([hN, hC, gN (1) gC (1)), {“平均控制组”，“卵巢癌组平均”，.．.“控制信封”，“卵巢癌包膜”}，.．.“位置”，“西北”)标题(“组平均数及组信封”）
观察一下，显然没有一个单一的特征可以完美地区分这两个群体。
主要特性
寻找重要特征的一个简单方法是假设每个M/Z值是独立的，并计算一个双向t检验。rankfeatures返回最显著的M/Z值的索引，例如按检验统计量的绝对值排序的100个索引。这种特征选择方法也被称为过滤方法，其中学习算法不涉及如何选择特征。
[feat,stat] = rankfeatures(Y,grp，“标准”，的tt，“数量”, 100);
的第一个输出rankfeatures可用于提取显著特征的M/Z值。
sig_mass = MZ(壮举);sig_Masses (1:7)%显示前七个
Ans = 1.0e+03 * 8.1009 8.1016 8.1024 8.1001 8.1032 7.7366 7.7359
的第二个输出rankfeatures是具有检验统计量绝对值的向量。你可以用yyaxis．
图;yyaxis左plot(MZ， [mean_N mean_C]);Ylim ([-1,20]) xlim([7950,8300]) title(“重要M/Z值”)包含(xAxisLabel);ylabel (yAxisLabel);yyaxis正确的情节(MZ,统计);ylim ([22]) ylabel (检验统计量的）;传奇({“平均控制组”，“卵巢癌组平均”，“测试数据”})
注意在高M/Z值但低强度(~8100 Da.)有显著区域。度量类可分性的其他方法在rankfeatures，如基于熵的Bhattacharyya，或经验接受者工作特征(ROC)曲线下的面积。
基于线性判别分析的盲分类
现在您已经确定了一些重要的特征，您可以使用这些信息对癌症和正常样本进行分类。由于样本数量较少，您可以使用20%坚持值运行交叉验证，以更好地估计分类器性能。cvpartition允许为不同类型的系统评估方法设置训练和测试指标，如hold-out、K-fold和Leave-M-Out。
Per_eval = 0.20;%的训练规模进行交叉验证rng (“默认”）;初始化随机生成器到相同的状态%用于生成已发布的示例CV = cvpartition(grp，“坚持”per_eval)
cv =保留交叉验证分区NumObservations: 216 NumTestSets: 1 TrainSize: 173 TestSize: 43
观察到特征只从训练子集中选择，并且验证是用测试子集执行的。classperf允许您跟踪多个验证。
Cp_lda1 = classperf(grp);初始化CP对象为k = 1:10运行交叉验证10次CV =重分区(CV);feat = rankfeatures(Y(:，training(cv))，grp(training(cv))，“数量”, 100);c =分类(Y(功绩,测试(cv)), Y(成绩、培训(cv)), grp(培训(cv)));classperf (cp_lda1 c测试(cv));使用当前验证更新CP对象结束
在循环之后，您可以使用CP对象中的任何属性来评估整体盲分类的性能，例如错误率、灵敏度、特异性等。
cp_lda1
Label: " Description: " ClassLabels: {2x1 cell} GroundTruth: [216x1 double] NumberOfObservations: 216 ControlClasses: 2 TargetClasses: 1 ValidationCounter: 10 SampleDistribution: [216x1 double] ErrorDistribution: [216x1 double] SampleDistributionByClass: [2x1 double] ErrorDistributionByClass: [2x1 double] CountingMatrix: [3x2 double] CorrectRate: 0.8488 ErrorRate: 0.1512 LastCorrectRate: 0.8837 LastErrorRate: 0.1163 InconclusiveRate: 0 ClassifiedRate: 1 Sensitivity: 0.8208 Specificity:0.8842 positivepredictivue: 0.8995 negativepredictivue: 0.7962 positivelikeue: 7.0890 negativelikeue: 0.2026患病率:0.5581诊断表:[2x2 double]
这种朴素的特征选择方法可以根据区域信息消除一些特征。例如，'NWEIGHT' inrankfeatures大于相邻M/Z特征的测试统计量，以便其他重要的M/Z值可以合并到所选特征的子集中
Cp_lda2 = classperf(grp);初始化CP对象为k = 1:10运行交叉验证10次CV =重分区(CV);feat = rankfeatures(Y(:，training(cv))，grp(training(cv))，“数量”, 100,“NWEIGHT”5);c =分类(Y(功绩,测试(cv)), Y(成绩、培训(cv)), grp(培训(cv)));classperf (cp_lda2 c测试(cv));使用当前验证更新CP对象结束cp_lda2。CorrectRate平均正确分类率
Ans = 0.9023
数据维数的PCA/LDA降维
Lilien等人在[2]中提出了一种使用主成分分析(PCA)来降低数据维数的算法，然后使用LDA对组进行分类。在本例中，M/Z空间中2000个最重要的特征被映射到150个主成分
Cp_pcalda = classperf(grp);初始化CP对象为k = 1:10运行交叉验证10次CV =重分区(CV);选择2000个最重要的特征。feat = rankfeatures(Y(:，training(cv))，grp(training(cv))，“数量”, 2000);% PCA来降低维数P = pca(Y(feat,training(cv))');投影到PCA空间x = Y(feat，:)' * P(:，1:150);%使用LDAc =分类(x(测试(cv):), x(培训(cv):), grp(培训(cv)));classperf (cp_pcalda c测试(cv));结束cp_pcalda。CorrectRate平均正确分类率
Ans = 0.9814
子集特征选择的随机搜索
特征选择也可以通过分类来加强，这种方法通常被称为包装器选择方法。随机搜索特征选择生成随机特征子集，并使用学习算法独立评估其质量。然后，它从最常见的好的特征中选择一个池。Li等人在[3]中将这一概念应用于蛋白质表达模式的分析。的randfeatures函数允许您在随机化的特征子集上使用LDA或k-最近邻分类器搜索特征子集。
注意:以下示例是计算密集型的，因此已从示例中禁用。此外，为了获得更好的结果，您应该从中的默认值增加分类器的池大小和严格性randfeatures．类型帮助randfeatures获取更多信息。
如果0% <==更改为1以启用。这可能需要很长时间才能完成。CV =重分区(CV);[feat,fCount] = randfeatures(Y(:，training(cv))，grp(training(cv))，.．.“分类”，“哒”，“PerformanceThreshold”, 0.90);其他的负载randFeatCancerDetect结束
用评估集评估所选特征的质量
的第一个输出randfeatures是MZ值的索引的有序列表。第一项在分类效果好的子集中出现频率最高。第二个输出是每个值被选择的实际次数。你可以使用嘘看看这个分布。
图;嘘(fCount,马克斯(fCount) + 1);
您将看到大多数值在选定的子集中最多出现一次。放大可以更好地了解更频繁选择的值的细节。
轴([0 80 0 100])
只有少数值被选择超过10次。您可以通过使用主干图来显示最常选择的特征来可视化这些特征。
图;持有在；sigfeat = fCount;sigfeat (sigfeat <=10) = 0;情节(MZ (mean_N mean_C]);茎(MZ (sigFeats > 0) sigFeats (sigFeats > 0),“r”）;轴([2000、12000、80])传说({“平均控制组”，“卵巢癌组平均”，的显著特征}，.．.“位置”，“西北”)包含(xAxisLabel);ylabel (yAxisLabel);
这些特征似乎聚集在几组中。您可以通过下面的实验进一步调查有多少特征是重要的。使用最频繁选择的特征对数据进行分类，然后使用两个最频繁选择的特征，依此类推，直到使用所有被选择超过10次的特征。然后，您可以看到添加更多特征是否会改善分类器。
nSig = sum(fCount>10);cp_rndfeat = 0 (20,nSig);为i = 1:nSig为J = 1:20 CV =重新划分(CV);P = pca(Y(feat(1:i)，training(cv))');x = Y(feat(1:i)，:)' * P;c =分类(x(测试(cv):), x(培训(cv):), grp(培训(cv)));Cp = classperf(grp,c,test(cv));cp_rndfeat(j,i) = cp.CorrectRate;平均正确分类率结束结束图(1:nSig， [max(cp_rndfeat);mean(cp_rndfeat)]);传奇({“最佳CorrectRate”，“意思是CorrectRate”}，“位置”，“东南”）
从这张图中你可以看到，即使只有三个特征，有时也可以得到完美的分类。你还会注意到，平均正确率的最大值出现在少数特征上，然后逐渐下降。
[bestAverageCR, bestNumFeatures] = max(mean(cp_rndfeat));
现在，您可以可视化给出最佳平均分类的特征。您可以看到，这些实际上只对应于数据中的三个峰值。
图;持有在；sigfeat = fCount;sigfeat (sigfeat <=10) = 0;ax_handle = plot(MZ，[mean_N mean_C]);茎(MZ(壮举(1:bestNumFeatures)), sigFeats(壮举(1:bestNumFeatures)),“r”）;轴([7650、8850、80])传说({“平均控制组”，“卵巢癌组平均”，的显著特征})包含(xAxisLabel);ylabel (yAxisLabel);
替代统计学习算法
MATLAB®中有许多分类工具，您也可以使用它们来分析蛋白质组学数据。其中包括支持向量机(fitcsvm)、k-最近邻(fitcknn)、神经网络(深度学习工具箱™)、分类树(fitctree)．对于特征选择，还可以使用顺序子集特征选择(sequentialfs)或使用遗传算法(全局优化工具箱)优化随机搜索方法。例如，请参见质谱数据特征的遗传算法搜索．
参考文献
[1]Conrads, T P, V A Fusaro, S Ross, D Johann, V Rajapakse, B A Hitt, S M Steinberg等，“用于卵巢癌检测的高分辨率血清蛋白质组学特征”。内分泌相关癌症，2004,6,163-78。
[２]莉莲，瑞安·H，哈尼·法里德，布鲁斯·r·唐纳德。人类血清质谱分析中表达依赖蛋白组数据的概率疾病分类计算生物学杂志第10期，no。6(2003年12月):925-46。
[3]李，L, D. M.乌姆巴赫，P.特里和J. A.泰勒。“GA/KNN方法在SELDI蛋白质组学数据中的应用”生物信息学，没有。10 (July 1,2004): 1638-40。
[4]Petricoin, Emanuel F, Ali M Ardekani, Ben A Hitt, Peter J Levine, Vincent A Fusaro, Seth M Steinberg, Gordon B Mills等，“利用血清蛋白质组学模式识别卵巢癌”。《柳叶刀》359，没有。9306(2002年2月):572-77。
另请参阅
相关的话题
您也可以从以下列表中选择一个网站: