В этом примере показаны рабочие процессы для анализа данных экспрессии генов с помощью алгоритма метагена-аттрактора. Данные по экспрессии генов доступны для многих модельных организмов и болезненных состояний. В этом примере показано, как использовать metafeatures функция для изучения биомолекулярных фенотипов при раке молочной железы.
Атлас генома рака (TCGA) включает несколько видов данных по нескольким показаниям рака. TCGA включает измерения экспрессии генов, экспрессии белков, клинических результатов и многое другое. В этом примере вы исследуете экспрессию генов рака молочной железы.
Исследователи собрали образцы опухолей и использовали микрочипы Agilent G4502A для измерения экспрессии их генов. В этом примере используются Level-3 данные выражения, которые были обработаны из исходных измерений в вызовы выражения. Данные были получены 20 мая 2014 года.
Загрузите данные в MATLAB ®. Файл MATTCGA_Breast_Gene_Expression.mat содержит данные экспрессии генов 17814 генов для 590 различных пациентов. Данные выражения хранятся в переменной geneExpression. Имена генов хранятся в переменной geneNames.
load TCGA_Breast_Gene_ExpressionЧтобы убедиться в организации данных, проверьте количество генов и образцов в этом наборе данных.
size(geneExpression)
ans = 1×2
17814 590
geneNames - клеточный массив имен генов. Доступ к записям можно получить с помощью индексирования массива ячеек MATLAB:
geneNames{655}ans = 'EGFR'
Этот массив ячеек указывает, что 655-я строка переменной geneExpression содержит измерения экспрессии для экспрессии гена рецептора эпидермального фактора роста (EGFR).
Алгоритм метагена-аттрактора был разработан в рамках задачи DREAM 8 по разработке прогностических биомаркеров для выживания при раке молочной железы. Подход метагена-аттрактора обнаруживает и количественно оценивает лежащие в основе биомолекулярные события. Эти события снижают размерность данных экспрессии генов, а также позволяют классифицировать подтипы и исследовать регуляторные механизмы [1].
Метаген определяется как любая взвешенная сумма экспрессии гена. Предположим, что у вас есть коллекция генов совместной экспрессии. Метаген можно создать путем усреднения уровней экспрессии генов в коллекции.
Существует потенциал для улучшения нашего понимания экспрессии генов, захваченных в этом метагене. Предположим, вы создадите набор весов, которые количественно определяют сходство между генами в нашей коллекции и метагеном. Гены, которые более похожи на метаген, получают больший вес, в то время как гены, которые менее похожи, получают меньший вес. Используя эти новые веса, можно сформировать новый метаген, который является средневзвешенным уровнем экспрессии генов. Новый метаген лучше фиксирует биомолекулярное событие, которое управляет некоторым элементом регуляции генов в данных экспрессии.
Эта процедура формирует ядро алгоритма метагена аттрактора. Сформируйте метаген, используя некоторую текущую оценку весов, затем обновите весовые коэффициенты на основе меры подобия. Метагены-аттракторы определяются как притягивающие фиксированные точки этого итеративного процесса.
Алгоритм существует в широком семействе неподконтрольных алгоритмов машинного обучения. Связанные алгоритмы включают анализ основных компонентов, различные алгоритмы кластеризации (особенно нечеткие c-средства), факторизацию неотрицательной матрицы и другие. Основным преимуществом метагенного подхода является то, что результаты алгоритма, как правило, более четко связаны с фенотипом, определяемым экспрессией генов.
Конкретно, в i-ой итерации алгоритма. У вас есть вектор весов, , размером 1 на число генов. Оценка метагена во время i-й итерации составляет:
* G
- количество генов по количеству образцов матрицы экспрессии генов. Для обновления весов:
Mi, Gj)
+ 1 - j-й элемент Wi 1, Gj - j-й ряд G, а J - метрика подобия. В алгоритме метагенного аттрактора J определяется как:
, Gj) α
если корреляция между и больше 0. MI - это взаимная информация между и . Функция metafeatures использует B-сплайн-оценщик взаимной информации, описанной в [3].
Если вместо этого корреляция между и меньше или равна 0, то:
= 0
Все веса больше или равны нулю. Поскольку взаимная информация инвариантна, вы можете нормализовать веса любым способом. Здесь они нормализованы, поэтому их сумма равна 1.
Алгоритм инициализируется случайными или выбранными пользователем весами. Он продолжается до тех пор, пока изменение между итерациями не станет небольшим или не будет исчерпано заранее определенное число итераций.
Данные имеют несколько NaN значения. Чтобы проверить, сколько, суммируйте показатель, возвращенный isnan.
sum(sum(isnan(geneExpression)))
ans = 1695
Из примерно 10 миллионов записей geneExpression, отсутствуют 1695 записей. Прежде чем продолжить, вам нужно будет разобраться с этими отсутствующими записями.
Существует несколько способов расчета этих отсутствующих значений. Можно использовать простой метод, называемый K ближайшего соседнего вменения, предоставленный с помощью панели инструментов биоинформатики (TM). K-ближайшее вменение соседа работает путем замены отсутствующих данных соответствующим значением из средневзвешенного значения k ближайшие столбцы к столбцу с отсутствующими данными.
Использовать k = 3 и замените текущее значение geneExpression с тем, у которого нет NaN значения.
geneExpression = knnimpute(geneExpression,3);
Переменная geneExpression не имеет NaN значения.
sum(sum(isnan(geneExpression)))
ans = 0
Для получения дополнительной информации о knnimpute, см. документацию по панели инструментов биоинформатики.
doc knnimpute
Функция metafeatures использует алгоритм метагена-аттрактора для идентификации мотивов регуляции генов.
Настройка структуры опций. В этом случае установите отображение для предоставления информации об алгоритме в каждой итерации.
opts = struct('Display','iter');
metafeatures также позволяет задавать начальные значения. Вы можете затравить стартовые веса, чтобы подчеркнуть гены, которые вам интересны. Существует три распространенных фактора рака молочной железы, ERBB2 (также называемый HER2), эстроген и прогестрон.
Установите вес для каждого из этих генов равным 1 в трех различных строках startValues. Каждая строка соответствует начальным значениям для различных копий. strcmp сравнивает интересующие гены и список генов в наборе данных. find возвращает индекс в списке гена.
erbb = find(strcmp('ERBB2',geneNames)); estrogen = find(strcmp('ESR1',geneNames)); progestrone = find(strcmp('PGR',geneNames)); startValues = zeros(size(geneExpression,1),3); startValues(erbb,1) = 1; startValues(estrogen,2) = 1; startValues(progestrone,3) = 1;
Звонить metafeatures с набором вмененных данных. Второй аргумент, geneNames - список всех генов в наборе данных. Ввод имен генов не требуется. Тем не менее, имена генов могут позволить исследовать высоко ранжированные гены, которые возвращаются алгоритмом, чтобы получить представление о биомолекулярном событии, описанном метагеном.
[meta, weights, genes_sorted] = metafeatures(geneExpression,geneNames,'start',startValues,'options',opts);
Caching self information ...
... done. Took 43.8561 seconds.
Caching entropy and binning information...
... done. Took 20.4185 seconds.
non-zero
Found iter diff weights
1 1 1.26e+01 8924
1 2 7.29e+00 8885
1 3 4.22e+00 8796
1 4 2.54e+00 8761
1 5 1.63e+00 8745
1 6 1.14e+00 8720
1 7 8.59e-01 8706
1 8 7.18e-01 8682
1 9 7.04e-01 8687
1 10 6.44e-01 8680
1 11 5.53e-01 8676
1 12 4.56e-01 8664
1 13 3.67e-01 8654
1 14 2.91e-01 8649
1 15 2.30e-01 8642
1 16 1.83e-01 8636
1 17 1.46e-01 8634
1 18 1.17e-01 8631
1 19 9.45e-02 8632
1 20 7.65e-02 8634
1 21 6.22e-02 8633
1 22 5.06e-02 8631
1 23 4.13e-02 8635
1 24 3.38e-02 8639
1 25 2.76e-02 8636
1 26 2.26e-02 8633
1 27 1.85e-02 8633
1 28 1.51e-02 8635
1 29 1.24e-02 8635
1 30 1.02e-02 8634
1 31 8.35e-03 8633
1 32 6.85e-03 8633
1 33 5.57e-03 8633
1 34 4.59e-03 8631
1 35 3.78e-03 8631
1 36 3.07e-03 8632
1 37 2.53e-03 8632
1 38 2.06e-03 8632
1 39 1.70e-03 8632
1 40 1.40e-03 8632
1 41 1.15e-03 8632
1 42 9.24e-04 8632
1 43 7.70e-04 8632
1 44 6.21e-04 8632
1 45 5.20e-04 8632
1 46 4.43e-04 8632
1 47 3.49e-04 8632
1 48 2.97e-04 8632
1 49 2.36e-04 8632
1 50 1.93e-04 8632
1 51 1.56e-04 8632
1 52 1.42e-04 8632
1 53 8.98e-05 8632
1 54 9.72e-05 8632
1 55 5.37e-05 8632
1 56 7.47e-05 8632
1 57 5.17e-05 8632
1 58 4.81e-05 8632
1 59 2.85e-05 8632
1 60 1.97e-05 8632
1 61 3.05e-05 8632
1 62 1.41e-05 8632
1 63 1.02e-05 8632
1 64 7.89e-06 8632
1 65 9.34e-06 8632
1 66 2.07e-05 8632
1 67 1.52e-05 8632
1 68 2.26e-05 8632
1 69 1.55e-05 8632
1 70 2.24e-05 8632
1 71 1.75e-05 8632
1 72 2.01e-05 8632
1 73 6.47e-06 8632
1 74 1.62e-05 8632
1 75 2.23e-05 8632
1 76 1.93e-05 8632
1 77 1.71e-05 8632
1 78 6.94e-06 8632
1 79 3.21e-06 8632
1 80 1.58e-05 8632
1 81 2.02e-05 8632
1 82 1.99e-05 8632
1 83 2.12e-05 8632
1 84 1.79e-05 8632
1 85 1.60e-05 8632
1 86 1.78e-05 8632
1 87 1.87e-05 8632
1 88 1.66e-05 8632
1 89 5.98e-06 8632
1 90 1.26e-05 8632
1 91 2.14e-05 8632
1 92 1.82e-05 8632
1 93 6.97e-06 8632
1 94 1.04e-05 8632
1 95 2.13e-05 8632
1 96 6.39e-06 8632
1 97 1.75e-05 8632
1 98 2.37e-05 8632
1 99 2.01e-05 8632
1 100 1.98e-05 8632
Warning: 'Maximum iterations exceeded, terminating early.'
2 1 1.93e+01 9893
2 2 6.04e+00 9885
2 3 3.80e+00 9883
2 4 2.53e+00 9886
2 5 1.73e+00 9881
2 6 1.13e+00 9873
2 7 7.19e-01 9869
2 8 4.63e-01 9866
2 9 3.08e-01 9870
2 10 2.13e-01 9874
2 11 1.54e-01 9872
2 12 1.15e-01 9874
2 13 8.72e-02 9874
2 14 6.68e-02 9874
2 15 5.14e-02 9874
2 16 3.97e-02 9875
2 17 3.07e-02 9875
2 18 2.37e-02 9873
2 19 1.84e-02 9871
2 20 1.42e-02 9871
2 21 1.10e-02 9871
2 22 8.54e-03 9872
2 23 6.62e-03 9872
2 24 5.05e-03 9872
2 25 4.01e-03 9872
2 26 3.09e-03 9872
2 27 2.38e-03 9872
2 28 1.85e-03 9872
2 29 1.43e-03 9872
2 30 1.09e-03 9872
2 31 8.46e-04 9872
2 32 6.73e-04 9872
2 33 5.10e-04 9872
2 34 3.81e-04 9872
2 35 2.98e-04 9872
2 36 2.46e-04 9872
2 37 1.51e-04 9872
2 38 1.63e-04 9872
2 39 1.15e-04 9872
2 40 7.11e-05 9872
2 41 1.18e-04 9872
2 42 7.28e-05 9872
2 43 1.89e-05 9872
2 44 4.24e-05 9872
2 45 1.60e-05 9872
2 46 6.75e-06 9872
2 47 4.81e-05 9872
2 48 2.47e-05 9872
2 49 1.04e-05 9872
2 50 7.46e-06 9872
2 51 9.31e-06 9872
2 52 5.25e-06 9872
2 53 3.89e-05 9872
2 54 9.38e-06 9872
2 55 3.33e-05 9872
2 56 1.48e-05 9872
2 57 2.45e-05 9872
2 58 2.58e-05 9872
2 59 1.00e-05 9872
2 60 1.86e-05 9872
2 61 5.87e-05 9872
2 62 2.97e-05 9872
2 63 1.07e-05 9872
2 64 8.84e-06 9872
2 65 8.29e-06 9872
2 66 1.58e-05 9872
2 67 1.48e-05 9872
2 68 5.00e-06 9872
2 69 2.74e-05 9872
2 70 1.20e-05 9872
2 71 2.91e-05 9872
2 72 9.45e-06 9872
2 73 1.75e-05 9872
2 74 1.56e-05 9872
2 75 6.56e-06 9872
2 76 1.79e-05 9872
2 77 2.67e-05 9872
2 78 5.55e-05 9872
2 79 2.55e-05 9872
2 80 1.03e-05 9872
2 81 2.74e-05 9872
2 82 2.04e-05 9872
2 83 1.00e-05 9872
2 84 1.11e-05 9872
2 85 9.83e-06 9872
2 86 2.71e-05 9872
2 87 1.42e-05 9872
2 88 1.28e-05 9872
2 89 2.24e-05 9872
2 90 4.58e-05 9872
2 91 3.36e-05 9872
2 92 9.74e-06 9872
2 93 1.06e-05 9872
2 94 1.50e-05 9872
2 95 5.05e-05 9872
2 96 1.12e-05 9872
2 97 2.52e-05 9872
2 98 9.77e-06 9872
2 99 6.10e-06 9872
2 100 2.97e-05 9872
Warning: 'Maximum iterations exceeded, terminating early.'
3 1 3.75e+00 9963
3 2 1.08e+00 9966
3 3 4.29e-01 9959
3 4 1.87e-01 9961
3 5 8.45e-02 9958
3 6 3.88e-02 9957
3 7 1.80e-02 9956
3 8 8.36e-03 9956
3 9 3.89e-03 9956
3 10 1.78e-03 9956
3 11 8.68e-04 9956
3 12 3.96e-04 9956
3 13 1.89e-04 9956
3 14 8.92e-05 9956
3 15 4.25e-05 9956
3 16 1.16e-05 9956
3 17 1.57e-05 9956
3 18 1.67e-05 9956
3 19 1.59e-05 9956
3 20 1.07e-05 9956
3 21 9.21e-06 9956
3 22 1.59e-05 9956
3 23 6.23e-06 9956
3 24 8.68e-06 9956
3 25 1.56e-05 9956
3 26 1.55e-05 9956
3 27 9.65e-06 9956
3 28 9.74e-06 9956
3 29 9.75e-06 9956
3 30 9.75e-06 9956
3 31 9.84e-06 9956
3 32 1.49e-05 9956
3 33 1.05e-05 9956
3 34 1.43e-05 9956
3 35 2.14e-05 9956
3 36 6.64e-06 9956
3 37 1.54e-06 9956
3 38 2.23e-06 9956
3 39 2.98e-06 9956
3 40 8.89e-06 9956
3 41 1.60e-05 9956
3 42 1.06e-05 9956
3 43 9.08e-06 9956
3 44 1.60e-05 9956
3 45 6.74e-06 9956
3 46 8.71e-06 9956
3 47 9.45e-06 9956
3 48 1.48e-05 9956
3 49 1.05e-05 9956
3 50 1.43e-05 9956
3 51 2.15e-05 9956
3 52 6.64e-06 9956
3 53 1.41e-06 9956
3 54 1.98e-06 9956
3 55 2.58e-06 9956
3 56 9.07e-06 9956
3 57 6.54e-06 9956
3 58 5.44e-06 9956
3 59 4.36e-06 9956
3 60 8.43e-06 9956
3 61 1.08e-05 9956
3 62 9.73e-06 9956
3 63 9.72e-06 9956
3 64 9.72e-06 9956
3 65 9.75e-06 9956
3 66 9.78e-06 9956
3 67 9.82e-06 9956
3 68 1.50e-05 9956
3 69 1.02e-05 9956
3 70 1.34e-05 9956
3 71 2.08e-05 9956
3 72 1.30e-05 9956
3 73 2.11e-05 9956
3 74 1.56e-05 9956
3 75 9.45e-06 9956
3 76 1.48e-05 9956
3 77 1.11e-05 9956
3 78 8.97e-06 9956
3 79 1.31e-05 9956
3 80 2.19e-05 9956
3 81 8.99e-06 9956
3 82 1.60e-05 9956
3 83 7.51e-06 9956
3 84 6.78e-06 9956
3 85 7.51e-06 9956
3 86 1.10e-05 9956
3 87 1.39e-05 9956
3 88 6.38e-06 9956
3 89 6.05e-06 9956
3 90 4.66e-06 9956
3 91 7.28e-06 9956
3 92 7.98e-06 9956
3 93 1.15e-05 9956
3 94 8.72e-06 9956
3 95 1.56e-05 9956
3 96 1.82e-05 9956
3 97 1.23e-05 9956
3 98 6.69e-06 9956
3 99 1.63e-06 9956
3 100 1.15e-06 9956
Warning: 'Maximum iterations exceeded, terminating early.'
Переменная meta имеет значение трех метагенов, обнаруженных для каждого образца. Вы можете построить графики трех метагенов, чтобы получить представление о природе регуляции генов в различных фенотипах рака молочной железы.
plot3(meta(1,:),meta(2,:),meta(3,:),'o') xlabel('ERBB2 metagene') ylabel('Estrogen metagene') zlabel('Progestrone metagene')
На сюжете можно наблюдать несколько вещей.
На сюжете присутствует группа точек, скомпонованных вместе с низкими значениями для всех трёх метагенов. Исходя из уровней мРНК, ожидается, что точки связаны с образцами опухоли, которые являются тройными отрицательными или базальными.
Существует также группа точек, которые имеют высокую экспрессию метагена рецептора эстрогена. Эта группа охватывает как высокую, так и низкую экспрессию метагена прогестрона. Нет точек с высокой экспрессией метагена прогестрона и низкой экспрессией метагена эстрогена. Это открытие согласуется с наблюдением, что рак молочной железы ER-/PR + встречается крайне редко [2].
Остальные моменты - это ERBB2 положительные виды рака. Они имеют меньшее представление в этом наборе данных, чем гормональный и тройной негативный рак. Также не существует прочных связей между экспрессией гормональных рецепторов и статусом ERBB2.
Чтобы лучше понять регуляцию генов, захваченную метагенами, более внимательно посмотрите на метаген, обнаруженный путем инициализации рецептора эстрогена с массой 1. Вы можете перечислить десять лучших генов, вносящих вклад в метаген для 11-го обнаруженного метагена.
genes_sorted(1:10,2)
ans = 10x1 cell
{'AGR3' }
{'ESR1' }
{'CA12' }
{'AGR2' }
{'MLPH' }
{'FOXA1'}
{'THSD4'}
{'FSIP1'}
{'ANXA9'}
{'XBP1' }
Этот метаген фиксирует биомолекулярное событие, связанное с переходом к эстрогеновому раку молочной железы. Перечисленные четыре высших гена:
Гомолог переднего градиента 3 (AGR3)
Рецептор эстрогена 1 (ESR1)
Угольная ангидраза 12 (CA12)
Гомолог переднего градиента 2 (AGR2)
Изменения транскрипции в каждом из этих генов вовлечены в рак молочной железы, вызванный эстрогеном. Известно, что три гена, отличные от ESR1, совместно экспрессируются с ESR1. Идентификация этих генов иллюстрирует способность алгоритма метагена-аттрактора связывать экспрессию генов с фенотипами.
Аналогичные варианты метагена эстрогена и метагена ERBB2 описаны в [1]. Порядок вкладов генов несколько отличается между этим анализом и [1], потому что использовался другой набор данных о раке молочной железы. Следует ожидать вариаций весов, но порядок генов по весам примерно одинаков. В частности, гены с топ-10 весами в основном одинаковы между этой версией, и версией, описанной в [1]. Аналогично, существует значительное перекрытие между генами с топ-100 весов.
Гены могут способствовать множеству метагенов. В этом смысле алгоритм метагена-аттрактора является «мягкой» техникой кластеризации. В этом примере, обнаруживая метагены в данных о раке молочной железы, наблюдается совпадение в наборах генов, которые имеют больший вес вклада в метагены эстрогена и прогестрона.
Если вес «повышен», когда он превышает 0,001, то:
elevated_weights = weights>.001;
Сумма столбцов elevated_weights - общее число повышенных весов в каждом из трех метагенов.
sum(elevated_weights)
ans = 1×3
19 96 27
Из 96 повышенных весов для метагена эстрогена и 27 для метагена прогестрона имеется 22 повышенных веса, которые находятся в обоих наборах.
sum(elevated_weights(:,2) & elevated_weights(:,3))
ans = 22
Однако между метагеном ERBB2 и метагеном эстрогена нет никакого перекрытия:
sum(elevated_weights(:,1) & elevated_weights(:,2))
ans = 0
а также отсутствие перекрытия между метагеном ERBB2 и метагеном прогестрона:
sum(elevated_weights(:,1) & elevated_weights(:,3))
ans = 0
В метрике подобия алгоритма параметр alpha управляет степенью нелинейности. По мере увеличения альфа количество метагенов имеет тенденцию к увеличению. Альфа-значение по умолчанию равно 5, поскольку это значение было хорошим для работы в [1], но для различных наборов данных или вариантов использования необходимо скорректировать альфа-значение.
Чтобы проиллюстрировать эффекты альфа, если альфа равен 1 в анализе рака молочной железы, то прогестерон и метагены эстрогена не отличаются.
[meta_alpha_1, weights_alpha_1, genes_sorted_alpha_1] = ... metafeatures(geneExpression,geneNames,'start',startValues,'alpha',1);
Warning: 'Maximum iterations exceeded, terminating early.'
Warning: 'Maximum iterations exceeded, terminating early.'
Warning: 'Maximum iterations exceeded, terminating early.'
В этом случае возвращаются только два метагена, несмотря на то, что мы трижды запускали алгоритм.
size(meta_alpha_1)
ans = 1×2
2 590
Этот результат объясняется тем, что по умолчанию metafeatures возвращает только уникальные метагены. Инициализация с весом для ESR1 установлена в 1, и инициализация с весом для PGR установлена в 1, оба сходятся к метагенам, которые фактически одинаковы.
[1] Чэн, Вэй-И, Тай-Хсиен Оу Ян и Димитрис Анастасиу. «Биомолекулярные события при раке, выявленные метагенами-аттракторами». Вычислительная биология PLoS 9.2 (2013): e1002920.
[2] Hefti, Marco M., et al. «Рецептор эстрогена отрицательный/рецептор прогестерона положительный рак молочной железы не является воспроизводимым подтипом». Исследование рака молочной железы 15.4 (2013): R68.
[3] Daub, Carsten O., et al. «Оценка взаимной информации с использованием функций B-сплайна? улучшенная мера сходства для анализа данных экспрессии генов». БМК биоинформатика 5.1 (2004): 118.