Идентификация биомолекулярных подгрупп с использованием Metagenes аттрактора
Этот пример показывает рабочие процессы для анализа данных экспрессии генов с помощью алгоритма аттрактора метагена. Данные по экспрессии генов доступны для многих модельных организмов и условий. В этом примере показано, как использовать metafeatures функция для исследования биомолекулярных фенотипов при раке молочной железы.
Данные атласа генома рака
Атлас генома рака (TCGA) включает несколько видов данных по нескольким признакам рака. TCGA включает измерения экспрессии генов, экспрессии белка, клинических результатов и многое другое. В этом примере вы исследуете экспрессию гена рака молочной железы.
Исследователи собрали выборки опухоли и использовали микромассивы Agilent G4502A для измерения экспрессии их генов. В этом примере вы используете данные выражения Level-3, которые были постобработаны из исходных измерений в вызовы выражения. Данные были получены 20 мая 2014 года.
Загрузите данные в MATLAB ®. Файл MAT TCGA_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 по количеству генов. Оценка метагена во время Вторая итерация:
количество генов по количеству выборок матрицы экспрессии генов. Для обновления весов:
является 1й элемент , является вторая строка G, и является метрикой подобия. В алгоритме аттрактора метагена, определяется как:
если корреляция между и больше 0. MI - это взаимная информация между и . Функция metafeatures использует B-сплайн оценку взаимной информации, описанную в [3].
Если вместо этого корреляция между и меньше или равно 0, тогда:
Веса все больше или равны нулю. Поскольку взаимная информация является инвариантной шкалы, можно нормализовать веса любым способом, который вы выбираете. Здесь они нормированы, поэтому их сумма равна 1.
Алгоритм инициализируется случайными или выбранными пользователем весами. Он продолжается до изменения в между итерациями мало, или предварительно заданное количество итераций исчерпано.
Очистка данных
Данные имеют несколько NaN значения. Чтобы проверить, сколько, сумма по индикатору, возвращенному isnan.
sum(sum(isnan(geneExpression)))
ans = 1695
Из примерно 10 миллионов записей geneExpression, отсутствуют 1695 записей. Прежде чем продолжить, вам нужно будет разобраться с этими отсутствующими записями.
Существует несколько способов вписать эти отсутствующие значения. Можно использовать простой метод, называемый K ближайшим соседом, предоставленным Bioinformatics Toolbox(TM). K-ближайших соседей вменение работает путем замены недостающих данных на соответствующее значение из взвешенного среднего значения k ближайшие столбцы к столбцу с отсутствующими данными.
Использование k = 3 и замените текущее значение geneExpression с тем, у кого нет NaN значения.
geneExpression = knnimpute(geneExpression,3);
Переменная geneExpression не имеет NaN значения.
sum(sum(isnan(geneExpression)))
ans = 0
Для получения дополнительной информации о knnimpute, см. документацию Bioinformatics Toolbox.
doc knnimpute
Идентификация биомолекулярных событий с помощью алгоритма метагена аттрактора
Функция metafeatures использует алгоритм аттрактора метагена, чтобы идентифицировать мотивы регуляции генов.
Setup структуры опций. В этом случае установите отображение, чтобы предоставить информацию об алгоритме при каждой итерации.
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 весами.
Гены могут способствовать нескольким метагенам. В этом смысле алгоритм метагена аттрактора является «мягким» методом кластеризации. В этом примере, находя метагены в данных о раке молочной железы, происходит перекрытие наборов генов, которые имеют большие веса вклада в метагены эстрогена и прогестрона.
Если вес «повышен», когда он больше .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
Роль Альфа
В метрике подобия алгоритма параметр альфа управляет степенью нелинейности. Когда альфа увеличивается, количество метагенов, как правило, увеличивается. Альфа по умолчанию является 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] Cheng, Wei-Yi, Tai-Hsien Ou Yang, and Dimitris Anastassiou. «Биомолекулярные события при раке, выявленные метагенами-аттракторами». Вычислительная биология PLoS 9.2 (2013): e1002920.
[2] Hefti, Marco M., et al. «Отрицательный рецептор эстрогена/рецептор прогестерона положительный рак молочной железы не является воспроизводимым подтипом». Исследование рака молочной железы 15.4 (2013): R68.
[3] Daub, Carsten O., et al. «Оценка взаимной информации с использованием B-сплайн функций? улучшенная мера подобия для анализа данных экспрессии генов». BMC биоинформатика 5.1 (2004): 118.