Введение

С платформами олигонуклеотида микромассивов 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

Фоновый метод корректировки 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

Процедура 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.