Визуализация и предварительная обработка написанных через дефис наборов данных масс-спектрометрии для профилирования метаболита и белка/Пептида

Этот пример показывает, как управлять, предварительно обработать и визуализировать данные из Жидкостной хроматографии вместе с Масс-спектрометрией (LC/MS). Эти большие и высокие размерные наборы данных экстенсивно используются в протеомике и metabolomics исследовании. Визуализация комплексного пептида или смесей метаболита предоставляет интуитивный метод, чтобы оценить демонстрационное качество. Кроме того, методическое исправление и предварительная обработка могут привести к автоматизированному высокому анализу пропускной способности выборок, позволяющих точную идентификацию значительных метаболитов и определенных функций пептида в биологической выборке.

Введение

В типичном написанном через дефис эксперименте масс-спектрометрии белки и метаболиты получены от ячеек, тканей, или биологических жидкостей, расторгли и денатурировали в решении, и ферментативным образом переварили в смеси. Эти смеси затем разделяются или Высокоэффективной жидкостной хроматографией (HPLC), капиллярным электрофорезом (CE) или газовой хроматографией (GC) и связываются с методом идентификации масс-спектрометрии, таким как Масс-спектрометрия Ионизации Электроспрея (ESI-MS), матрица помогла ионизации (MALDI или TOF-MS SELDI), или тандемная масс-спектрометрия (MS/MS).

Открытые Репозитории данных и mzXML Формат файла

В данном примере мы используем тестовую выборку набор данных LC-ESI-MS с семью соединениями белка. Данные в этом примере (7MIX_STD_110802_1) от Репозитория данных сашими. Набор данных не распределяется с MATLAB®. Чтобы завершить этот пример, необходимо загрузить набор данных в локальную директорию или собственный репозиторий. Также можно попробовать другие наборы данных, доступные в других общедоступных базах данных для данных о выражении белка, таких как Атлас Пептида в Институте Системной биологии [1].

Большинство текущих массовых спектрометров может перевести или сохранить данные о приобретении с помощью mzXML схемы. Этот стандарт является XML (расширяемый язык разметки) - базирующийся формат общего файла, разработанный проектом сашими обратиться к проблемам, вовлеченным в представление наборов данных от различных производителей и от различных экспериментальных настроек в общую и расширяемую схему. файлы mzXML, используемые в написанной через дефис масс-спектрометрии, являются обычно очень большими. Функция MZXMLINFO позволяет вам получать основную информацию о наборе данных, не читая его в память. Например, можно получить количество сканирований, область значений времени задержания и количество тандемных инструментов MS (уровни).

info = mzxmlinfo('7MIX_STD_110802_1.mzXML','NUMOFLEVELS',true)
info = 

  struct with fields:

                            Filename: '7MIX_STD_110802_1.mzXML'
                         FileModDate: '01-Feb-2013 11:54:30'
                            FileSize: 26789612
                       NumberOfScans: 7161
                           StartTime: 'PT0.00683333S'
                             EndTime: 'PT200.036S'
       DataProcessingIntensityCutoff: 'N/A'
            DataProcessingCentroided: 'false'
            DataProcessingDeisotoped: 'N/A'
    DataProcessingChargeDeconvoluted: 'N/A'
       DataProcessingSpotIntegration: 'N/A'
                    NumberOfMSLevels: 2

Функция MZXMLREAD читает XML-документ в структуру MATLAB. Поля scan и index помещаются в первый уровень выходной структуры для улучшенного доступа к спектральным данным. Остаток от mzXML дерева документов анализируется согласно спецификациям схемы. Этот набор данных LC/MS содержит 7 161 сканирование с двумя уровнями MS. Для этого примера вы будете использовать только первые сканирования уровня. Вторые спектры уровня обычно используются для идентификации пептида/белка и прибывают в более поздний этап в некоторых типах исследований рабочего процесса. MZXMLREAD может отфильтровать желаемые сканирования, не загружая весь набор данных в память:

mzXML_struct = mzxmlread('7MIX_STD_110802_1.mzXML','LEVEL',1)
Starting to parse document...
Building mzXML substructure...
Building scan substructure...
Building index substructure...
DONE!

mzXML_struct = 

  struct with fields:

     scan: [2387x1 struct]
    mzXML: [1x1 struct]
    index: [1x1 struct]

Если вы получаете какие-либо ошибки, связанные с памятью или пространством "кучи" Java во время загрузки, попытайтесь увеличить свое пространство "кучи" Java, как описано здесь.

Более подробная информация, принадлежащая массовый спектрометр и экспериментальные условия, найдена в поле msRun.

mzXML_struct.mzXML.msRun
ans = 

  struct with fields:

         scanCount: 7161
         startTime: 'PT0.00683333S'
           endTime: 'PT200.036S'
        parentFile: [1x1 struct]
      msInstrument: [1x1 struct]
    dataProcessing: [1x1 struct]

Чтобы упростить обработку данных, функция MZXML2PEAKS извлекает список peaks от каждого сканирования в массив ячеек (peaks]) и их соответствующее время задержания в вектор-столбец (time). Можно извлечь спектры определенного уровня путем установки параметра входа LEVEL.

[peaks,time] = mzxml2peaks(mzXML_struct);
numScans = numel(peaks)
numScans =

        2387

Функция MSDOTPLOT создает отображение обзора самого интенсивного peaks в целом наборе данных. В этом случае мы визуализируем только самые интенсивные 5% ионного peaks интенсивности путем установки входного параметра QUANTILE на 0,95.

h = msdotplot(peaks,time,'quantile',.95);
title('5 Percent Overall Most Intense Peaks')

Можно также отфильтровать peaks индивидуально для каждого сканирования с помощью процентили основной пиковой интенсивности. Основной пик является самым интенсивным пиком, найденным на каждом сканировании [2]. Этот параметр дан автоматически большинством спектрометров. Эта операция требует, чтобы запрос в mxXML структуру получил основную пиковую информацию. Обратите внимание на то, что вы могли также найти основную пиковую интенсивность путем итерации функции MAX по пиковому списку.

basePeakInt = [mzXML_struct.scan.basePeakIntensity]';
peaks_fil = cell(numScans,1);
for i = 1:numScans
    h = peaks{i}(:,2) > (basePeakInt(i).*0.75);
    peaks_fil{i} = peaks{i}(h,:);
end

whos('basePeakInt','level_1','peaks','peaks_fil')
msdotplot(peaks_fil,time)
title('Peaks Above (0.75 x Base Peak Intensity) for Each Scan')
  Name                Size               Bytes  Class     Attributes

  basePeakInt      2387x1                19096  double              
  peaks            2387x1             14050896  cell                
  peaks_fil        2387x1               308664  cell                

Можно настроить 3-D обзор отфильтрованного peaks с помощью функции STEM3. Функция STEM3 требует, чтобы поместить данные в три вектора, элементы которых формируют триплеты (время задержания, масса/заряд и значение интенсивности), которые представляют каждую основу.

peaks_3D = cell(numScans,1);
for i = 1:numScans
    peaks_3D{i}(:,[2 3]) = peaks_fil{i};
    peaks_3D{i}(:,1) = time(i);
end
peaks_3D = cell2mat(peaks_3D);

figure
stem3(peaks_3D(:,1),peaks_3D(:,2),peaks_3D(:,3),'marker','none')
axis([0 12000 400 1500 0 1e9])
view(60,60)
xlabel('Retention Time (seconds)')
ylabel('Mass/Charge (M/Z)')
zlabel('Relative Ion Intensity')
title('Peaks Above (0.75 x Base Peak Intensity) for Each Scan')

Можно построить окрашенные основы с помощью функции PATCH. Для каждого триплета в peaks_3D чередуйте новый триплет с набором значений интенсивности, чтобы обнулить. Затем создайте цветной вектор, зависящий от интенсивности основы. Логарифмическое преобразование улучшает динамический диапазон палитры. Поскольку чередованные триплеты присваивают NaN, так, чтобы функция PATCH не проводила линии, соединяющие непрерывные основы.

peaks_patch = sortrows(repmat(peaks_3D,2,1));
peaks_patch(2:2:end,3) = 0;

col_vec = log(peaks_patch(:,3));
col_vec(2:2:end) = NaN;

figure
patch(peaks_patch(:,1),peaks_patch(:,2),peaks_patch(:,3),col_vec,...
      'edgeColor','flat','markeredgecolor','flat','Marker','x','MarkerSize',3);
axis([0 12000 400 1500 0 1e9])
view(90,90)
xlabel('Retention Time (seconds)')
ylabel('Mass/Charge (M/Z)')
zlabel('Relative Ion Intensity')
title('Peaks Above (0.75 x Base Peak Intensity) for Each Scan')

view(40,40)

Создание Карт Тепла Наборов данных LC/MS

Общие методы в промышленности работают с пиковой информацией (иначе центроидные данные) вместо необработанных сигналов. Это может сохранить память, но некоторые важные детали не видимы, особенно когда необходимо осмотреть выборки с комплексными смесями. Чтобы далее анализировать этот набор данных, мы можем создать общую сетку в размерности массы/заряда. С тех пор не все сканирования имеют достаточно информации, чтобы восстановить исходный сигнал, мы используем пиковый метод передискретизации сохранения. Путем выбора соответствующих параметров для функции MSPPRESAMPLE можно гарантировать, что разрешение спектров не потеряно, и что максимальные значения peaks коррелируют к исходной пиковой информации.

[MZ,Y] = msppresample(peaks,5000);
whos('MZ','Y')
  Name         Size                 Bytes  Class     Attributes

  MZ        5000x1                  20000  single              
  Y         5000x2387            47740000  single              

С этой матрицей ионной интенсивности, Y, можно создать цветную карту тепла. Функция MSHEATMAP автоматически настраивает шкалу палитры, используемую, чтобы показать статистически значительный peaks с горячими цветами и шумный peaks с холодными цветами. Алгоритм основан на кластеризирующемся значительном peaks и шумном peaks путем оценки смеси Gaussians с подходом Максимизации Ожидания. Кроме того, можно использовать параметр входа MIDPOINT, чтобы выбрать произвольный порог, чтобы разделить шумный peaks от значительного peaks, или можно в интерактивном режиме переключить палитру, чтобы скрыть или вывести на экран peaks. При работе с картами тепла распространено отобразить логарифм ионной интенсивности, который улучшает динамический диапазон палитры.

fh1 = msheatmap(MZ,time,log(Y),'resolution',.1,'range',[500 900]);
title('Original LC/MS Data Set')

Можно масштабировать к маленьким видимым областям, чтобы наблюдать данные, или в интерактивном режиме или программно использование функции AXIS. Мы наблюдаем некоторые области с высокой относительной ионной интенсивностью. Они представляют пептиды в биологической выборке.

axis([527 539 385 496])

Можно наложить исходную пиковую информацию набора данных LC/MS. Это позволяет вам оценить производительность операции передискретизации сохранения пика. Можно использовать возвращенный указатель, чтобы скрыться/вывести на экран точки.

dp1 = msdotplot(peaks,time);

Калибровка Местоположения Массы/Заряда Peaks к Общей Сетке

Двумерный peaks, кажется, является шумным, или они не показывают компактную форму в непрерывных спектрах. Это - типичная проблема для многих массовых спектрометров. Случайные колебания значения массы/заряда, извлеченного от peaks, реплицируют профили, как, сообщают, колеблются от 0,1% до 0,3% [3]. Такая изменчивость может быть вызвана несколькими факторами, например, плохой калибровкой детектора, низкого отношения сигнал-шум или проблем в пиковых алгоритмах экстракции. Функция MSPALIGN реализует усовершенствованные алгоритмы раскладывания данных, которые синхронизируют все спектры в наборе данных к общей сетке массы/заряда (CMZ). CMZ может быть выбран произвольно, или он может быть оценен после анализа данных [2,4,5]. Пиковая процедура соответствия может использовать или самое близкое соседнее правило или динамическое выравнивание программирования.

[CMZ, peaks_CMZ] = mspalign(peaks);

Повторите процесс визуализации с выровненным peaks: выполните пиковую передискретизацию сохранения, создайте карту тепла, наложите выровненную пиковую информацию и изменение масштаба в ту же видимую область как прежде. Когда спектр перекалиброван, возможно отличить изотопные шаблоны некоторых пептидов.

[MZ_A,Y_A] = msppresample(peaks_CMZ,5000);
fh2 = msheatmap(MZ_A,time,log(Y_A),'resolution',.10,'range',[500 900]);
title('LC/MS Data Set with the Mass/Charge Calibrated to a CMZ')
dp2 = msdotplot(peaks_CMZ,time);
axis([527 539 385 496])

Калибровка Местоположения Массы/Заряда Peaks Локально

MSPALIGN вычисляет один CMZ для целого набора данных LC/MS. Это не может быть идеальным случаем для выборок с более комплексными смесями пептидов и/или метаболитов, чем набор данных, используемый в этом примере. В случае комплексных смесей можно выровнять каждый спектр к локальному набору спектров, которые содержат только информативный peaks (высокая интенсивность) с подобными временами задержания, в противном случае калибровочный процесс в областях с маленьким peaks (низкая интенсивность) может быть смещен другим peaks, который совместно использует подобные значения массы/заряда, но является в различные времена задержания. Чтобы выполнить более прекрасную калибровку, можно использовать функцию SAMPLEALIGN. Эта функция является обобщением алгоритмов Ограниченного динамического времени, деформируясь (CDTW), обычно используемых в речи, обрабатывающей [6]. В следующем цикле for мы поддерживаем буфер с интенсивностью предыдущих выровненных спектров (LAI). Ионная интенсивность спектров масштабируется с анонимной функцией SFSAMPLEALIGN), чтобы уменьшить расстояние между большим peaks. SAMPLEALIGN уменьшает общее расстояние всех совпадающих точек и вводит разрывы по мере необходимости. В этом случае мы используем более прекрасный вектор MZ (FMZ), такой, что мы сохраняем правильное значение массы/заряда peaks как можно больше. Примечание: это может занять время, когда алгоритм CDTW выполняется 2,387 раз.

SF = @(x) 1-exp(-x./5e7); % scaling function
DF = @(R,S) sqrt((SF(R(:,2))-SF(S(:,2))).^2 + (R(:,1)-S(:,1)).^2);

FMZ = (500:0.15:900)';  % setup a finer MZ vector
LAI = zeros(size(FMZ)); % init buffer for the last alignment intensities

peaks_FMZ = cell(numScans,1);
for i = 1:numScans
    % show progress
    if ~rem(i,250)
        fprintf(' %d...',i);
    end
    % align peaks in current scan to LAI
    [k,j] = samplealign([FMZ,LAI],double(peaks{i}),'band',1.5,'gap',[0,2],'dist',DF);
    % updating the LAI buffer
    LAI = LAI*.25;
    LAI(k) = LAI(k) + peaks{i}(j,2);
    % save the alignment
    peaks_FMZ{i} = [FMZ(k) peaks{i}(j,2)];
end
 250... 500... 750... 1000... 1250... 1500... 1750... 2000... 2250...

Повторите процесс визуализации и изменение масштаба к видимой области.

[MZ_B,Y_B] = msppresample(peaks_FMZ,4000);
fh3 = msheatmap(MZ_B,time,log(Y_B),'resolution',.10,'range',[500 900]);
title('LC/MS Data Set with the Mass/Charge Calibrated Locally')
dp3 = msdotplot(peaks_FMZ,time);
axis([527 539 385 496])

Как последний шаг, чтобы улучшить изображение, можно применить Гауссов фильтр в хроматографическом направлении, чтобы сглаживать целый набор данных.

Gpulse = exp(-.1*(-10:10).^2)./sum(exp(-.1*(-10:10).^2));
YF = convn(Y_B,Gpulse,'same');
fh4 = msheatmap(MZ_B,time,log(YF),'resolution',.10,'limits',[500 900]);
title('Final Enhanced LC/MS Data Set')
dp4 = msdotplot(peaks_FMZ,time);
axis([527 539 385 496])

Можно соединить оси двух карт тепла, к в интерактивном режиме или программно сравнить области между двумя наборами данных. В этом случае мы сравниваем оригинал, и финал улучшил матрицы LC/MS.

linkaxes(findobj([fh1 fh4],'Tag','MSHeatMap'))
axis([521 538 266 617])

Извлечение спектров Используя общую ионную хроматограмму

Если набор данных LC/MS сглаживается и передискретизируется в обычную сетку, возможно извлечь самые информативные спектры путем рассмотрения локальных максимумов Общей ионной хроматограммы (TIC). TIC прямо вычисляется путем подведения итогов строк YF. Затем используйте функцию MSPEAKS, чтобы найти временные стоимости задержания для извлечения выбранных подмножеств спектров.

TIC = mean(YF);
pt = mspeaks(time,TIC','multiplier',10,'overseg',100,'showplot',true);
title('Extracted peaks from the Total Ion Chromatogram (TIC)')
pt(pt(:,1)>4000,:) = []; % remove spectra above 4000 seconds
numPeaks = size(pt,1)
numPeaks =

    22

Создайте 3-D график выбранных спектров.

xRows = samplealign(time,pt(:,1),'width',1); % finds the time index for every peak
xSpec = YF(:,xRows);                         % gets the signals to plot

figure;
hold on
box on
plot3(repmat(MZ_B,1,numPeaks),repmat(1:numPeaks,numel(MZ_B),1),xSpec)
view(20,85)

ax = gca;
ax.YTick = 1:numPeaks;
ax.YTickLabel = num2str(time(xRows));
axis([500 900 0 numPeaks 0 1e8])

xlabel('Mass/Charge (M/Z)')
ylabel('Time')
zlabel('Relative Ion Intensity')
title('Extracted Spectra Subset')

Наложите маркеры для извлеченных спектров по расширенной тепловой карте.

linkaxes(findobj(fh4,'Tag','MSHeatMap'),'off')
figure(fh4)
hold on
for i = 1:numPeaks
    plot([400 1500],xRows([i i]),'m')
end
axis([500 900 100 925])
dp4.Visible = 'off';
title('Final Enhanced LC/MS Data Set with Extracted Spectra Marked')

Ссылки

[1] Desiere, F. и др., "Проект Атласа Пептида", Исследование Нуклеиновых кислот, 34:D655-8, 2006.

[2] Purvine, S., Kolker, N. и Kolker, E., "Спектральная качественная оценка для тандемной протеомики масс-спектрометрии Высокой Пропускной способности", OMICS: журнал интегральной биологии, 8 (3):255-65, 2004.

[3] Kazmi, A.S., и др., "Выравнивание спектров массы высокого разрешения: Разработка эвристического подхода для metabolomics", Metabolomics, 2 (2):75-83, 2006.

[4] Jeffries, N., "Алгоритмы для выравнивания масс-спектрометрии протеомные данные", Биоинформатика, 21 (14):3066-3073, 2005.

[5] Ю, W., и др., "Несколько достигают максимума выравнивание в последовательном анализе данных: основанный на пробеле шкалы подход", Сделка IEEE®/ACM Вычислительная Биология и Биоинформатика, 3 (3):208-219, 2006.

[6] Sakoe, H. и Чиба s., "Динамическая оптимизация алгоритма программирования для распознавания произносимого слова", Сделка IEEE. Акустика, Речь и Обработка сигналов, ASSP-26 (1):43-9, 1978.