Этот пример показывает, как использовать MATLAB® и Bioinformatics Toolbox™ для предварительной обработки данных тестового уровня олигонуклеотида Affymetrix® микромассивов с двумя методами предварительной обработки, Устойчивым средним значением мультимассивов (RMA) и GC Устойчивое Среднее значение Мультимассивов (GCRMA).
С платформами олигонуклеотида микромассивов Affymetrix экспрессия гена измеряется с помощью тестовых наборов, состоящих из 11 - 20 идеальных пар (PM), зонды (25 нуклеотидов в длине) дополнительный, чтобы предназначаться для mRNA последовательностей. Каждый тестовый набор также имеет то же количество несоответствия (MM) зонды, в которых 13-й нуклеотид был изменен на свое дополнение. Зонды Премьер-министра разработаны для гена определенная гибридизация. Измерения зонда MM управления, как думают, включают большую часть фона неопределенная привязка, такая как перекрестная гибридизация. Зонд ПРЕМЬЕР-МИНИСТРА и его соответствующий зонд MM упоминаются как тестовая пара.
Измеренная тестовая интенсивность и местоположения от гибридизированного микромассива хранятся в файле CEL. Для каждой платформы Affymetrix микромассивов информация, связывающая тестовые пары, чтобы зондировать идентификаторы набора, и к местоположениям на массиве, хранится в файле библиотеки CDF. Тестовая информация о последовательности хранится в файле последовательности (FASTA или разделенный от вкладки формат).
В целом предварительная обработка данных о выражении тестового уровня Affymetrix состоит из трех шагов: фоновая корректировка, нормализация и резюмирование на тестовом уровне набора как мера уровня экспрессии соответствующего mRNA. Много методов существуют для статистических процедур этих трех шагов. Два популярных метода, RMA (Irizarry и др., 2003) и GCRMA (Ву и др., 2004), используются в этом примере.
Примечание: Этот пример показывает RMA и GCRMA предварительную обработку процедур, чтобы вычислить значения выражения из файлов входа CEL в постепенных деталях, с помощью нескольких функций. Можно также завершить тот же RMA или методы GCRMA в одном вызове функции при помощи Bioinformatics Toolbox affyrma или функции affygcrma, соответственно.
Общедоступный набор данных, содержащий измерения Affymetrix микромассивов 42 тканей опухоли эмбриональной центральной нервной системы (CNS, Pomeroy и др., 2002), используется для этого примера. Вы импортируете и получите доступ к тестовым данным об уровне нескольких массивов, и затем выполните измерения уровня экспрессии с RMA и GCRMA предварительная обработка методов.
Эксперимент CNS проводился с помощью массива Affymetrix HuGeneFL GeneChip®, и данные хранились в файлах CEL. Информация, связанная с каждым зондом, содержится в файле библиотеки Affymetrix Hu6800 CDF.
Если у вас уже нет файла библиотеки Hu6800 CDF, загрузите zip-файл библиотеки HuGeneFL Genome Array. Извлеките Hu6800. CDF-файл в директорию, такую как C:\Examples\affypreprocessdemo\libfiles. Примечание: необходимо будет указать в порядке получить доступ к файлам библиотеки, но вы не должны запускать файл setup.exe в архиве.
Набор данных CNS (файлы CEL) доступен здесь. Чтобы завершить этот пример, загрузите файлы CEL набора данных CNS в директорию, такие как C:\Examples\affypreprocessdemo\data. Разархивируйте архивы файла CEL. Примечание: Этот набор данных содержит больше файлов CEL, чем необходимо для этого примера.
CNS_DataA_Sample_CEL.txt, файл, которому предоставляют Bioinformatics Toolbox, содержит список 42 имен файлов CEL, используемых для этого примера и выборок (10 medulloblastomas, 10 rhabdoid, 10 злокачественных глиом, 8 супратенториальных PNETS и 4 нормальных человеческих мозжечка), которому они принадлежат. Загрузите эти данные в две переменные MATLAB.
fid = fopen('CNS_DataA_Sample_CEL.txt','r'); ftext = textscan(fid,'%q%q'); fclose(fid); samples = ftext{1}; cels = ftext{2};
Установите переменные celPath и libPath на пути файлов CEL и директорий библиотеки.
celPath = 'C:\Examples\affypreprocessdemo\data'; libPath = 'C:\Examples\affypreprocessdemo\libfiles';
Переименуйте файлы буфера перемещаемого изображения так, чтобы каждое имя файла запустилось с номера MG, который следует за подчеркиванием "_" в исходном имени файла. Например, GSM1688666_MG1999060202AA.CEL переименован к MG1999060202AA.CEL. Вы не должны запускать этот код, если имена файлов уже находятся в требуемом формате.
A = dir(fullfile(celPath,'*.cel')); fileNames = string({A.name}); for iFile = 1:numel(A) newName = fullfile(celPath,extractAfter(fileNames(iFile),"_")); movefile(fullfile(celPath,fileNames(iFile)),newName); end
Функциональный celintensityread может считать несколько файлов CEL и получить доступ к файлу библиотеки CDF. Это возвращает структуру MATLAB, содержащую тестовую информацию и тестовую интенсивность. Матрицы PM и интенсивности MM из нескольких файлов CEL хранятся в полях PMIntensities и MMIntensities. В каждой тестовой матрице интенсивности индексы столбца соответствуют порядку, в котором были считаны файлы CEL, и каждая строка соответствует зонду. Создайте структуру MATLAB PM и интенсивности зонда MM путем загрузки данных из файлов CEL из директории, где файлы CEL хранятся, и передача в пути туда, где вы хранили файл библиотеки CDF. (Примечание: celintensityread сообщит о продвижении к окну команды MATLAB. Можно выключить отчет о выполнении работ путем установки входного параметра VERBOSE на ложь.)
probeData = celintensityread(cels, 'Hu6800.CDF',... 'celpath', celPath, 'cdfpath', libPath, 'pmonly', false)
Reading CDF file: Hu6800.CDF
Reading file 1 of 42: MG2000040501AA
Reading file 2 of 42: MG2000040502AA
Reading file 3 of 42: MG2000040504AA
Reading file 4 of 42: MG2000040505AA
Reading file 5 of 42: MG2000040508AA
Reading file 6 of 42: MG2000040509AA
Reading file 7 of 42: MG2000040510AA
Reading file 8 of 42: MG2000040511AA
Reading file 9 of 42: MG2000040512AA
Reading file 10 of 42: MG2000040513AA
Reading file 11 of 42: MG2000051201AA
Reading file 12 of 42: MG2000051202AA
Reading file 13 of 42: MG2000051204AA
Reading file 14 of 42: MG2000051205AA
Reading file 15 of 42: MG2000051209AA
Reading file 16 of 42: MG2000071102AA
Reading file 17 of 42: MG2000051207AA
Reading file 18 of 42: MG2000051208AA
Reading file 19 of 42: MG2000051211AA
Reading file 20 of 42: MG2000051213AA
Reading file 21 of 42: MG2000061902AA
Reading file 22 of 42: MG2000061903AA
Reading file 23 of 42: MG2000061904AA
Reading file 24 of 42: MG2000061905AA
Reading file 25 of 42: MG2000061906AA
Reading file 26 of 42: MG2000070709AA
Reading file 27 of 42: MG2000070710AA
Reading file 28 of 42: MG2000070711AA
Reading file 29 of 42: MG2000070712AA
Reading file 30 of 42: MG2000070713AA
Reading file 31 of 42: MG1999112206AA
Reading file 32 of 42: MG2000033109AA
Reading file 33 of 42: MG2000033106AA
Reading file 34 of 42: MG2000033107AA
Reading file 35 of 42: MG1999112202AA
Reading file 36 of 42: MG1999112204AA
Reading file 37 of 42: MG2000011801AA
Reading file 38 of 42: MG2000031503AA
Reading file 39 of 42: MG2000032015AA
Reading file 40 of 42: MG2000030308AA
Reading file 41 of 42: MG2000011803AA
Reading file 42 of 42: MG2000011807AA
probeData =
struct with fields:
CDFName: 'Hu6800.CDF'
CELNames: {1×42 cell}
NumChips: 42
NumProbeSets: 7129
NumProbes: 140983
ProbeSetIDs: {7129×1 cell}
ProbeIndices: [140983×1 uint8]
GroupNumbers: [140983×1 uint8]
PMIntensities: [140983×42 single]
MMIntensities: [140983×42 single]
Определите количество загруженных файлов CEL.
nSamples = probeData.NumChips
nSamples =
42
Определите количество тестовых наборов на массиве HuGeneFL.
nProbeSets = probeData.NumProbeSets
nProbeSets =
7129
Определите количество зондов на массиве HuGeneFL.
nProbes = probeData.NumProbes
nProbes =
140983
Чтобы выполнить предварительную обработку GCRMA, тестовая информация последовательности массива HuGeneFL также запрашивается. Сайт поддержки Affymetrix предоставляет тестовую информацию о последовательности для большинства доступных массивов, или как FASTA отформатированные или файлы с разделением табуляцией. Этот пример принимает, что у вас есть файл HuGeneFL_probe_tab в директории файлов библиотеки. Используйте функциональный affyprobeseqread, чтобы проанализировать файл последовательности и возвратить тестовые последовательности в nProbes x 25 матриц целых чисел, которые представляют основы последовательности зонда премьер-министра со строками, соответствующими зондам на чипе и столбцах, соответствующих основным положениям 25-mer.
S = affyprobeseqread('HuGeneFL_probe_tab', 'Hu6800.CDF',... 'seqpath', libPath, 'cdfpath', libPath, 'seqonly', true)
S =
struct with fields:
SequenceMatrix: [140983×25 uint8]
Процедура RMA использует только интенсивность зонда PM для фоновой корректировки (Irizarry и др., 2003), в то время как GCRMA настраивает фон с помощью тестовой информации о последовательности и интенсивности зонда управления MM, чтобы оценить неопределенную привязку (Ву и др., 2004). И RMA и GCRMA предшествует нормализация квантиля (Bolstad и др., 2003) и среднее резюмирование полировки (Irizarry и др., 2003) интенсивности PM.
Фоновый метод корректировки RMA исправляет чип интенсивности зонда PM чипом. Интенсивность зонда Премьер-министра моделируется как сумма нормального шумового компонента и компонента экспоненциального сигнала. Использование rmabackadj к фону настраивает интенсивность премьер-министра в данных о CNS. Можно осмотреть гистограмму распределения интенсивности и предполагаемую фоновую корректировку определенного чипа путем установки входного параметра SHOWPLOT на индекс столбца чипа.
pms_bg = rmabackadj(probeData.PMIntensities, 'showplot', 1);
К нескольким нелинейным методам нормализации успешно применились микроданные массива Affymetrix. Процедура RMA нормирует данные тестового уровня с методом нормализации квантиля. Используйте quantilenorm, чтобы нормировать настроенную интенсивность PM фона в данных о CNS. Примечание: Если вы интересуетесь инвариантным рангом методом нормализации набора, используйте функцию affyinvarsetnorm вместо этого.
pms_bgnorm = quantilenorm(pms_bg);
Средняя процедура полировки применяется к интенсивности премьер-министра в резюмировании. Чтобы вычислить значения выражения, используйте rmasummary, чтобы обобщить тестовую интенсивность каждого тестового набора через несколько микросхем. Значения выражения являются тестовыми сводными данными интенсивности набора на журнале 2 шкалы.
cns_rma_exp = rmasummary(probeData.ProbeIndices, pms_bgnorm);
Процедура GCRMA настраивает для оптической шумовой и неопределенной привязки (NSB), учитывающей эффект более сильного связывания пар G/C (Naef и др., 2003, Ву и др., 2004). GCRMA использует тестовую информацию о последовательности, чтобы оценить тестовое сродство для вычисления неопределенной привязки. Тестовое сродство моделируется как сумма позиционно-зависимых основных эффектов. Обычно, тестовое сродство оценивается от интенсивности MM эксперимента NSB. Если данные NSB не доступны, тестовое сродство может все еще быть оценено от информации о последовательности и интенсивности зонда MM, нормированной тестовой интенсивностью медианы набора (Naef и др., 2003).
Для набора данных CNS используйте данные из микромассива, гибридизированного с нормальной выборкой мозжечков (Brain_Ncer_1), чтобы вычислить тестовое сродство для массива HuGeneFL. Используйте affyprobeaffinities, чтобы оценить тестовое сродство микромассива Affymetrix. Используйте параметр входа SHOWPLOT, чтобы осмотреть график, показывающий эффекты основы A, C, G, и T в этих 25 положениях.
figure idx = find(strcmpi('Brain_Ncer_1', samples)); [pmAlpha, mmAlpha] = affyprobeaffinities(S.SequenceMatrix,... probeData.MMIntensities(:, idx), 'showplot', true);
Примечание: существует 496 зондов на массиве HuGeneFL, которые не имеют информации о последовательности; сродством для этих зондов был NaN.
С тестовым доступным сродством сумма NSB может быть оценена путем подбора кривой кривой LOWESS через интенсивность зонда MM по сравнению со сродством зонда MM. Функциональный gcrmabackadj выполняет оптические и исправления NSB. Входной параметр, который SHOWPLOT показывает графику оптического шума, настроил интенсивность MM против своего сродства и сглаженный припадок заданного чипа. Можно вычислить фоновую интенсивность с одним из двух методов оценки, Оценки наибольшего правдоподобия (MLE) и Эмпирически-байесового (EB), который вычисляет следующее среднее значение определенной обязательной данной предшествующей наблюдаемой интенсивности. Здесь вы будете фон настраивать четыре массива с помощью обоих методов оценки. (Примечание: gcrmabackadj сообщит о продвижении к окну команды MATLAB. Можно выключить отчет о выполнении работ путем установки входного параметра VERBOSE на ложь.)
Фон настраивает первые четыре микросхем с помощью метода GCRMA-MLE и осматривает график интенсивности по сравнению со сродством к данным из третьего массива.
pms_MLE_bg = gcrmabackadj(probeData.PMIntensities(:,1:4),... probeData.MMIntensities(:, 1:4),... pmAlpha, mmAlpha, 'showplot', 3); % Adjust YLIM for better view ylim([4 16]);
Adjusting background for chip # 1 of 4 using MLE method. Adjusting background for chip # 2 of 4 using MLE method. Adjusting background for chip # 3 of 4 using MLE method. Adjusting background for chip # 4 of 4 using MLE method.
Фон настраивает первые четыре микросхем с помощью метода GCRMA-EB. Обработка с этим методом более в вычислительном отношении интенсивна и займет больше времени.
pms_EB_bg = gcrmabackadj(probeData.PMIntensities(:,1:4),... probeData.MMIntensities(:, 1:4),... pmAlpha, mmAlpha, 'method', 'EB');
Adjusting background for chip # 1 of 4 using EB method. Adjusting background for chip # 2 of 4 using EB method. Adjusting background for chip # 3 of 4 using EB method. Adjusting background for chip # 4 of 4 using EB method.
Можно продолжить предварительную обработку с quatilenorm и функциями rmasummary, или использовать функцию gcrma, чтобы сделать все. Функция gcrma выполняет фоновую корректировку и возвращается, меры по выражению фона настроили интенсивность зонда PM с помощью той же нормализации и методов резюмирования как RMA. Можно также передать в матрице последовательности вместо сродства. Функция автоматически вычислит сродство в этом случае. (Примечание: gcrma сообщит о продвижении к окну команды MATLAB. Можно выключить отчет о выполнении работ путем установки входного параметра VERBOSE на ложь.)
cns_mle_exp = gcrma(probeData.PMIntensities, probeData.MMIntensities,...
probeData.ProbeIndices, pmAlpha, mmAlpha);
Adjusting background for chip # 1 of 42 using MLE method. Adjusting background for chip # 2 of 42 using MLE method. Adjusting background for chip # 3 of 42 using MLE method. Adjusting background for chip # 4 of 42 using MLE method. Adjusting background for chip # 5 of 42 using MLE method. Adjusting background for chip # 6 of 42 using MLE method. Adjusting background for chip # 7 of 42 using MLE method. Adjusting background for chip # 8 of 42 using MLE method. Adjusting background for chip # 9 of 42 using MLE method. Adjusting background for chip # 10 of 42 using MLE method. Adjusting background for chip # 11 of 42 using MLE method. Adjusting background for chip # 12 of 42 using MLE method. Adjusting background for chip # 13 of 42 using MLE method. Adjusting background for chip # 14 of 42 using MLE method. Adjusting background for chip # 15 of 42 using MLE method. Adjusting background for chip # 16 of 42 using MLE method. Adjusting background for chip # 17 of 42 using MLE method. Adjusting background for chip # 18 of 42 using MLE method. Adjusting background for chip # 19 of 42 using MLE method. Adjusting background for chip # 20 of 42 using MLE method. Adjusting background for chip # 21 of 42 using MLE method. Adjusting background for chip # 22 of 42 using MLE method. Adjusting background for chip # 23 of 42 using MLE method. Adjusting background for chip # 24 of 42 using MLE method. Adjusting background for chip # 25 of 42 using MLE method. Adjusting background for chip # 26 of 42 using MLE method. Adjusting background for chip # 27 of 42 using MLE method. Adjusting background for chip # 28 of 42 using MLE method. Adjusting background for chip # 29 of 42 using MLE method. Adjusting background for chip # 30 of 42 using MLE method. Adjusting background for chip # 31 of 42 using MLE method. Adjusting background for chip # 32 of 42 using MLE method. Adjusting background for chip # 33 of 42 using MLE method. Adjusting background for chip # 34 of 42 using MLE method. Adjusting background for chip # 35 of 42 using MLE method. Adjusting background for chip # 36 of 42 using MLE method. Adjusting background for chip # 37 of 42 using MLE method. Adjusting background for chip # 38 of 42 using MLE method. Adjusting background for chip # 39 of 42 using MLE method. Adjusting background for chip # 40 of 42 using MLE method. Adjusting background for chip # 41 of 42 using MLE method. Adjusting background for chip # 42 of 42 using MLE method. Normalizing. Calculating expression.
Используйте коробчатые диаграммы, чтобы осмотреть дистрибутивы интенсивности премьер-министра первых четырех микросхем с тремя фоновыми процедурами корректировки.
figure subplot(4,1,1) maboxplot(log2(probeData.PMIntensities(:, 1:4)), samples(1:4),... 'title','Raw Intensities', 'orientation', 'horizontal') subplot(4,1,2) maboxplot(log2(pms_bg(:,1:4)), samples(1:4),... 'title','After RMA Background Adjustment','orient','horizontal') subplot(4,1,3) maboxplot(log2(pms_MLE_bg), samples(1:4),... 'title','After GCRMA-MLE Background Adjustment','orient','horizontal') subplot(4,1,4) maboxplot(log2(pms_EB_bg), samples(1:4),... 'title','After GCRMA-EB Background Adjustment','orient','horizontal')
Используйте коробчатые диаграммы, чтобы осмотреть фон исправленные и нормированные дистрибутивы интенсивности PM первых четырех микросхем с тремя фоновыми процедурами корректировки.
pms_MLE_bgnorm = quantilenorm(pms_MLE_bg); pms_EB_bgnorm = quantilenorm(pms_EB_bg); figure subplot(3,1,1) maboxplot(log2(pms_bgnorm(:, 1:4)), samples(1:4),... 'title','Normalized RMA Background Adjusted Intensity',... 'orientation', 'horizontal') subplot(3,1,2) maboxplot(log2(pms_MLE_bgnorm), samples(1:4),... 'title','Normalized GCRMA-MLE Background Adjusted Intensity',... 'orientation', 'horizontal') subplot(3,1,3) maboxplot(log2(pms_EB_bgnorm), samples(1:4),... 'title','Normalized GCRMA-EB Background Adjusted Intensity',... 'orientation', 'horizontal')
Можно выполнить импорт данных из файлов CEL и всех трех шагов предварительной обработки RMA и методов GCRMA, показанных в этом примере при помощи affyrma и функций affygcrma соответственно.
Для получения дополнительной информации об анализе экспрессии гена с Bioinformatics Toolbox смотрите Данные об Экспрессии гена Исследования Микромассивов.
Affymetrix и GeneChip являются зарегистрированными торговыми марками Affymetrix, Inc.
[1] Pomeroy, S.L., и др., "Прогноз центральной нервной системы эмбриональный результат опухоли на основе экспрессии гена", Природа, 415 (6870):436-42, 2002.
[2] Irizarry, R.A., и др., "Исследование, нормализация и сводные данные массива олигонуклеотида высокой плотности зондируют данные об уровне", Биостатистика, 4 (2):249-64, 2003.
[3] Ву, Z., и др., "Основанная на модели фоновая корректировка к массивам выражения олигонуклеотида", Журнал американской Статистической Ассоциации, 99 (468):909-17, 2004.
[4] Bolstad, B.M., и др., "Сравнение методов нормализации для данных массива олигонуклеотида высокой плотности на основе отклонения и смещения", Биоинформатика, 19 (2):185-93, 2003.
[5] Naef, F., и Magnasco, M.O. "Решая загадку ярких несоответствий: маркируя и эффективная привязка в массивах олигонуклеотида", Физический Анализ, E, Статистическая, Нелинейная и Мягкая Физика Вопроса, 68 (1Pt1):011906, 2003.