В этом примере показано, как управлять, предварительно обработайте и визуализируйте данные из Жидкостной хроматографии вместе с Масс-спектрометрией (LC/MS). Эти большие и высокие размерные наборы данных экстенсивно используются в протеомике и metabolomics исследовании. Визуализация комплексного пептида или смесей метаболита предоставляет интуитивный метод, чтобы оценить демонстрационное качество. Кроме того, методическая коррекция и предварительная обработка могут привести к автоматизированному высокому анализу пропускной способности выборок, позволяющих точную идентификацию значительных метаболитов и определенных функций пептида в биологической выборке.
В типичном написанном через дефис эксперименте масс-спектрометрии белки и метаболиты получены от ячеек, тканей, или биологических жидкостей, расторгли и денатурировали в решении, и ферментативным образом переварили в смеси. Эти смеси затем разделяются или Высокоэффективной жидкостной хроматографией (HPLC), капиллярным электрофорезом (CE) или газовой хроматографией (GC) и связываются с методом идентификации масс-спектрометрии, таким как Масс-спектрометрия Ионизации Электроспрея (ESI-MS), матрица помогла ионизации (MALDI или TOF-MS SELDI), или тандемная масс-спектрометрия (MS/MS).
В данном примере мы используем тестовую выборку набор данных 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
помещаются в первый уровень структуры output для улучшенного доступа к спектральным данным. Остаток от 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)
Общие методы в промышленности работают с пиковой информацией (иначе центроидные данные) вместо необработанных сигналов. Это может сохранить память, но некоторые важные детали не отображаются, особенно когда необходимо смотреть выборки с комплексными смесями. Чтобы далее анализировать этот набор данных, мы можем создать общую сетку в размерности массы/заряда. С тех пор не все сканирования имеют достаточно информации, чтобы восстановить исходный сигнал, мы используем пиковый метод передискретизации сохранения. Путем выбора соответствующих параметров для 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, реплицируют профили, как, сообщают, лежат в диапазоне от 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])
MSPALIGN
вычисляет один CMZ для целого набора данных LC/MS. Это не может быть идеальным случаем для выборок с более комплексными смесями пептидов и/или метаболитов, чем набор данных, используемый в этом примере. В случае комплексных смесей можно выровнять каждый спектр к локальному набору спектров, которые содержат только информативный peaks (высокая интенсивность) с подобными временами задержания, в противном случае калибровочный процесс в областях с маленьким peaks (низкая интенсивность) может быть смещен другим peaks, который совместно использует подобные значения массы/заряда, но является в различные времена задержания. Чтобы выполнить более прекрасную калибровку, можно использовать SAMPLEALIGN
функция. Эта функция является обобщением алгоритмов Ограниченной динамической трансформации временной шкалы (CDTW), обычно используемых в речи, обрабатывающей [6]. В следующем for
цикл, мы обеспечиваем буфер с интенсивностью предыдущих выровненных спектров (LAI). Ионная интенсивность спектров масштабируется с анонимной функцией SF
(в SAMPLEALIGN
) уменьшать расстояние между большим 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.