Работа с ГЕО серийными данными

Этот пример показывает, как получить ряд данных об экспрессии гена (GSE) из Автобуса экспрессии гена (GEO) NCBI и выполнить основной анализ профилей выражения.

Введение

Автобус экспрессии гена (GEO) NCBI является самым большим общедоступным репозиторием экспериментальных данных высокой пропускной способности микромассивов. ГЕО данные имеют четыре типа сущности включая ГЕО Платформу (GPL), ГЕО Выборка (GSM), ГЕО Ряд (GSE) и курировали ГЕО DataSet (GDS).

Запись Платформы описывает список элементов на массиве в эксперименте (например, комплементарные ДНК, олигонуклеотид проокружает). Каждая запись Платформы присвоена уникальный и стабильный ГЕО инвентарный номер (GPLxxx).

Демонстрационная запись описывает условия, при которых отдельная Выборка была обработана, манипуляции, которым это подверглось, и измерение распространенности каждого элемента, выведенного от него. Каждая Демонстрационная запись присвоена уникальный и стабильный ГЕО инвентарный номер (GSMxxx).

Серийная запись задает группу связанных Выборок и обеспечивает центр и описание целого исследования. Серийные записи могут также содержать таблицы, описывающие извлеченные данные, итоговые заключения или исследования. Каждая Серийная запись присвоена уникальный ГЕО инвентарный номер (GSExxx).

Запись DataSet (GDSxxx) представляет курировавший набор биологически и статистически сопоставимые ГЕО Выборки. ГЕО DataSets (GDSxxx) являются курировавшими наборами ГЕО Выборочных данных.

Больше информации о ГЕО может быть найдено в ГЕО Обзоре. Bioinformatics Toolbox™ обеспечивает функции, которые могут получить и проанализировать ГЕО файлы данных формата. GSE, GSM, GSD и данные GPL могут быть получены при помощи функции getgeodata. Функция getgeodata может также сохранить полученные данные в текстовом файле. ГЕО Серийные записи доступны в файлах формата SOFT и в файлах формата разделенного текста вкладки. Функциональный geoseriesread читает ГЕО Серийный файл текстового формата. Функция geosoftread читает обычно довольно большие файлы формата SOFT.

В этом примере вы получите набор данных GSE5847 из базы данных GEO и выполните статистическое тестирование на данных. ГЕО Серия GSE5847 содержит экспериментальные данные от исследования экспрессии гена основы опухоли и ячеек эпителия от 15 случаев воспалительного рака молочной железы (IBC) и 35 неподстрекательских случаев рака молочной железы (Boersma и др. 2008).

Получение ГЕО серийных данных

Функциональный getgeodata возвращает структуру, содержащую полученные данные от базы данных GEO. Можно также сохранить возвращенные данные в его исходном формате к локальной файловой системе для использования на последующих сеансах MATLAB®. Обратите внимание на то, что данные в общедоступных репозиториях часто курируются и обновляются; поэтому результаты этого примера могут немного отличаться, когда вы используете актуальные наборы данных.

gseData = getgeodata('GSE5847', 'ToFile', 'GSE5847.txt')

Используйте функцию geoseriesread, чтобы проанализировать ранее загруженный файл текстового формата GSE.

gseData = geoseriesread('GSE5847.txt')
gseData = 

  struct with fields:

    Header: [1x1 struct]
      Data: [22283x95 bioma.data.DataMatrix]

Возвращенная структура содержит поле Header, которое хранит метаданные Серийных данных и поле Data, которое хранит Серийные матричные данные.

Исследование данные GSE

Матричные данные GSE5847 в поле Data хранятся как объект DataMatrix. Объект DataMatrix является структурой данных как MATLAB 2D массив, но с дополнительными метаданными имен строки и имен столбцов. К свойствам DataMatrix можно получить доступ как другие объекты MATLAB.

get(gseData.Data)
            Name: ''
        RowNames: {22283x1 cell}
        ColNames: {1x95 cell}
           NRows: 22283
           NCols: 95
           NDims: 2
    ElementClass: 'double'

Имена строки являются идентификаторами тестовых наборов на массиве; имена столбцов являются ГЕО Демонстрационными инвентарными номерами.

gseData.Data(1:5, 1:5)
ans = 

                 GSM136326    GSM136327    GSM136328    GSM136329    GSM136330
    1007_s_at     10.45       9.3995       9.4248       9.4729       9.2788   
    1053_at      5.7195       4.8493       4.7321       4.7289       5.3264   
    117_at       5.9387       6.0833        6.448       6.1769       6.5446   
    121_at       8.0231       7.8947        8.345       8.1632       8.2338   
    1255_g_at    3.9548       3.9632       3.9641       4.0878       3.9989   

Серийные метаданные хранятся в поле Header. Структура содержит информацию о Ряде в поле Header.Series и демонстрационную информацию в поле Header.Sample.

gseData.Header
ans = 

  struct with fields:

     Series: [1x1 struct]
    Samples: [1x1 struct]

Поле ряда содержит заголовок эксперимента и микромассива ГЕО ID Платформы.

gseData.Header.Series
ans = 

  struct with fields:

                         title: 'Tumor and stroma from breast by LCM'
                 geo_accession: 'GSE5847'
                        status: 'Public on Sep 30 2007'
               submission_date: 'Sep 15 2006'
              last_update_date: 'Nov 14 2012'
                     pubmed_id: '17999412'
                       summary: 'Tumor epithelium and surrounding stromal cells were isolated using laser capture microdissection of human breast cancer to examine differences in gene expression based on tissue types from inflammatory and non-inflammatory breast cancer...'
                overall_design: 'We applied LCM to obtain samples enriched in tumor epithelium and stroma from 15 IBC and 35 non-IBC cases to study the relative contribution of each component to the IBC phenotype and to patient survival. '
                          type: 'Expression profiling by array'
                   contributor: 'Stefan,,Ambs...'
                     sample_id: 'GSM136326 GSM136327 GSM136328 GSM136329 GSM136330 GSM136331 GSM136332 GSM136333 GSM136334 GSM136335 GSM136336 GSM136337 GSM136338 GSM136339 GSM136340 GSM136341 GSM136342 GSM136343 GSM136344 GSM136345 GSM136346 GSM136347 GSM136348 GSM136349 GSM136350 GSM136351 GSM136352 GSM136353 GSM136354 GSM136355 GSM136356 GSM136357 GSM136358 GSM136359 GSM136360 GSM136361 GSM136362 GSM136363 GSM136364 GSM136365 GSM136366 GSM136367 GSM136368 GSM136369 GSM136370 GSM136371 GSM136372 GSM136373 GSM136374 GSM136375 GSM136376 GSM136377 GSM136378 GSM136379 GSM136380 GSM136381 GSM136382 GSM136383 GSM136384 GSM136385 GSM136386 GSM136387 GSM136388 GSM136389 GSM136390 GSM136391 GSM136392 GSM136393 GSM136394 GSM136395 GSM136396 GSM136397 GSM136398 GSM136399 GSM136400 GSM136401 GSM136402 GSM136403 GSM136404 GSM136405 GSM136406 GSM136407 GSM136408 GSM136409 GSM136410 GSM136411 GSM136412 GSM136413 GSM136414 GSM136415 GSM136416 GSM136417 GSM136418 GSM136419 GSM136420 '
                  contact_name: 'Stefan,,Ambs'
            contact_laboratory: 'LHC'
             contact_institute: 'NCI'
               contact_address: '37 Convent Dr Bldg 37 Room 3050'
                  contact_city: 'Bethesda'
                 contact_state: 'MD'
    contact_zip0x2Fpostal_code: '20892'
               contact_country: 'USA'
            supplementary_file: 'ftp://ftp.ncbi.nlm.nih.gov/pub/geo/DATA/supplementary/series/GSE5847/GSE5847_RAW.tar'
                   platform_id: 'GPL96'
                platform_taxid: '9606'
                  sample_taxid: '9606'
                      relation: 'BioProject: http://www.ncbi.nlm.nih.gov/bioproject/97251'

gseData.Header.Samples
ans = 

  struct with fields:

                         title: {1x95 cell}
                 geo_accession: {1x95 cell}
                        status: {1x95 cell}
               submission_date: {1x95 cell}
              last_update_date: {1x95 cell}
                          type: {1x95 cell}
                 channel_count: {1x95 cell}
               source_name_ch1: {1x95 cell}
                  organism_ch1: {1x95 cell}
           characteristics_ch1: {2x95 cell}
                  molecule_ch1: {1x95 cell}
          extract_protocol_ch1: {1x95 cell}
                     label_ch1: {1x95 cell}
            label_protocol_ch1: {1x95 cell}
                     taxid_ch1: {1x95 cell}
                  hyb_protocol: {1x95 cell}
                 scan_protocol: {1x95 cell}
                   description: {1x95 cell}
               data_processing: {1x95 cell}
                   platform_id: {1x95 cell}
                  contact_name: {1x95 cell}
            contact_laboratory: {1x95 cell}
             contact_institute: {1x95 cell}
               contact_address: {1x95 cell}
                  contact_city: {1x95 cell}
                 contact_state: {1x95 cell}
    contact_zip0x2Fpostal_code: {1x95 cell}
               contact_country: {1x95 cell}
            supplementary_file: {1x95 cell}
                data_row_count: {1x95 cell}

Поле data_processing содержит информацию методов предварительной обработки, в этом случае процедуры Устойчивого среднего значения мультимассивов (RMA).

gseData.Header.Samples.data_processing(1)
ans =

  1x1 cell array

    {'RMA'}

Поле source_name_ch1 содержит демонстрационный источник:

sampleSources = unique(gseData.Header.Samples.source_name_ch1);
sampleSources{:}
ans =

    'human breast cancer stroma'


ans =

    'human breast cancer tumor epithelium'

Поле Header.Samples.characteristics_ch1 содержит характеристики выборок.

gseData.Header.Samples.characteristics_ch1(:,1)
ans =

  2x1 cell array

    {'IBC'     }
    {'Deceased'}

Определите IBC и метки не-IBC для выборок от поля Header.Samples.characteristics_ch1 как метки группы.

sampleGrp = gseData.Header.Samples.characteristics_ch1(1,:);

Получение ГЕО платформы (GPL) данные

Серийные метаданные сказали нам ID платформы массивов: GPL96, который является массивом Affymetrix® GeneChip® Human Genome U133, установил HG-U133A. Получите файл формата GPL96 SOFT из ГЕО использования функции getgeodata. Например, примите, что вы использовали функцию getgeodata, чтобы получить файл Платформы GPL96 и сохраненный это к файлу, такому как GPL96.txt. Используйте функцию geosoftread, чтобы проанализировать этот файл формата SOFT.

gplData = geosoftread('GPL96.txt')
gplData = 

  struct with fields:

                 Scope: 'PLATFORM'
             Accession: 'GPL96'
                Header: [1x1 struct]
    ColumnDescriptions: {16x1 cell}
           ColumnNames: {16x1 cell}
                  Data: {22283x16 cell}

Поле ColumnNames структуры gplData содержит имена столбцов для данных:

gplData.ColumnNames
ans =

  16x1 cell array

    {'ID'                              }
    {'GB_ACC'                          }
    {'SPOT_ID'                         }
    {'Species Scientific Name'         }
    {'Annotation Date'                 }
    {'Sequence Type'                   }
    {'Sequence Source'                 }
    {'Target Description'              }
    {'Representative Public ID'        }
    {'Gene Title'                      }
    {'Gene Symbol'                     }
    {'ENTREZ_GENE_ID'                  }
    {'RefSeq Transcript ID'            }
    {'Gene Ontology Biological Process'}
    {'Gene Ontology Cellular Component'}
    {'Gene Ontology Molecular Function'}

Можно получить тестовые идентификаторы набора и названия генов для тестовых наборов платформы GPL69.

gplProbesetIDs = gplData.Data(:, strcmp(gplData.ColumnNames, 'ID'));
geneSymbols = gplData.Data(:, strcmp(gplData.ColumnNames, 'Gene Symbol'));

Используйте названия генов, чтобы маркировать гены в объекте DataMatrix gseData.Data. Следует иметь в виду, что набор зонда IDs из файла платформы не может быть в том же порядке как в gseData.Data. В этом примере они находятся в том же порядке.

Поменяйте имена строки данных о выражении к названиям генов.

gseData.Data = rownames(gseData.Data, ':', geneSymbols);

Отобразите первые пять строк и пять столбцов данных о выражении с именами строки как названия генов.

gseData.Data(1:5, 1:5)
ans = 

              GSM136326    GSM136327    GSM136328    GSM136329    GSM136330
    DDR1       10.45       9.3995       9.4248       9.4729       9.2788   
    RFC2      5.7195       4.8493       4.7321       4.7289       5.3264   
    HSPA6     5.9387       6.0833        6.448       6.1769       6.5446   
    PAX8      8.0231       7.8947        8.345       8.1632       8.2338   
    GUCA1A    3.9548       3.9632       3.9641       4.0878       3.9989   

Анализ данных

Bioinformatics Toolbox и Statistics and Machine Learning Toolbox™ предлагают широкий спектр инструментов анализа и визуализации для анализа микроданных массива. Однако, потому что это не наша главная цель объяснить методы анализа в этом примере, вы примените только несколько функций к профилю экспрессии гена от стромальных ячеек. Для более тщательно продуманных примеров об анализе экспрессии гена и доступных инструментах выбора функции, смотрите Данные об Экспрессии гена Исследования Микромассивов и Выбор Features for Classifying High-dimensional Data (Statistics and Machine Learning Toolbox).

Эксперимент выполнялся на IBC и выборках не-IBC, выведенных от стромальных ячеек и эпителиальных клеток. В этом примере вы будете работать с профилем экспрессии гена от стромальных ячеек. Определите демонстрационные индексы для стромального типа ячейки от поля gseData.Header.Samples.source_name_ch1.

stromaIdx = strcmpi(sampleSources{1}, gseData.Header.Samples.source_name_ch1);

Определите количество выборок от стромальных ячеек.

nStroma = sum(stromaIdx)
nStroma =

    47

stromaData = gseData.Data(:, stromaIdx);
stromaGrp = sampleGrp(stromaIdx);

Определите количество IBC и не-IBC стромальные выборки ячейки.

nStromaIBC = sum(strcmp('IBC', stromaGrp))
nStromaIBC =

    13

nStromaNonIBC = sum(strcmp('non-IBC', stromaGrp))
nStromaNonIBC =

    34

Можно также маркировать столбцы в stromaData с метками группы:

stromaData = colnames(stromaData, ':', stromaGrp);

Отобразите гистограмму нормированных измерений экспрессии гена заданного гена. Оси X представляют нормированный уровень экспрессии. Например, осмотрите распределение значений экспрессии гена этих генов.

fID = 331:339;
zValues = zscore(stromaData.(':')(':'), 0, 2);
bw = 0.25;
edges = -10:bw:10;
bins = edges(1:end-1) + diff(edges)/2;
histStroma = histc(zValues(fID, :)', edges) ./ (stromaData.NCols*bw);
figure;
for i = 1:length(fID)
    subplot(3,3,i);
    bar(edges, histStroma(:,i), 'histc')
    xlim([-3 3])
    if i <= length(fID)-3
        ax = gca;
        ax.XTickLabel = [];
    end
    title(sprintf('gene%d - %s', fID(i), stromaData.RowNames{fID(i)}))
end
suptitle('Gene Expression Value Distributions')

К профилю экспрессии гена получили доступ с помощью Affymetrix GeneChip более затем 22 000 функций на небольшом количестве выборок (~100). Среди 47 опухолей стромальные выборки существует 13 IBC и 34 не-IBC. Но не все гены дифференцированно выражаются между опухолями не-IBC и IBC. Статистические тесты необходимы, чтобы идентифицировать подпись экспрессии гена, которые отличают IBC от не-IBC стромальные выборки.

Используйте genevarfilter, чтобы отфильтровать гены с небольшим отклонением через выборки.

[mask, stromaData] = genevarfilter(stromaData);
stromaData.NRows
ans =

       20055

Примените t-статистическую-величину на каждый ген и сравните p-значения для каждого гена, чтобы найти значительно дифференцированно выраженные гены между IBC и группами не-IBC путем перестановки выборок (1,000 раз для этого примера).

rng default
[pvalues, tscores] = mattest(stromaData(:, 'IBC'), stromaData(:, 'non-IBC'),...
                           'Showhist', true', 'showplot', true, 'permute', 1000);

Выберите гены в заданном p-значении.

sum(pvalues < 0.001)
ans =

    52

Существует приблизительно 50 генов, выбранных непосредственно в p-значениях <0.001.

Сортировка и список лучшие 20 генов:

testResults = [pvalues, tscores];
testResults = sortrows(testResults);
testResults(1:20, :)
ans = 

                        p-values      t-scores
    INPP5E              2.3337e-05     5.0389 
    ARFRP1 /// IGLJ3    2.7575e-05     4.9753 
    USP46               3.4346e-05    -4.9054 
    GOLGB1              4.7706e-05    -4.7928 
    TTC3                0.00010702    -4.5053 
    THUMPD1             0.00013185    -4.4317 
                        0.00016066     4.3656 
    MAGED2              0.00017047    -4.3444 
    DNAJB9              0.00017833    -4.3266 
    KIF1C               0.00022129     4.2504 
                        0.00022251    -4.2482 
    DZIP3               0.00022425    -4.2454 
    COPB1               0.00023203    -4.2332 
    PSD3                0.00024659    -4.2138 
    PLEKHA4             0.00026506      4.186 
    DNAJB9              0.00027687    -4.1708 
    CNPY2               0.00028039    -4.1672 
    USP9X               0.00028454    -4.1619 
    SEC22B              0.00030187    -4.1392 
    GFER                0.00030527    -4.1352 

Ссылки

[1] Boersma, B.J., Reimers, M., И, M., Людвиг, J.A., и др. "Стромальная генная подпись сопоставила с воспалительным раком молочной железы", Международный журнал Рака, 122 (6):1324-32, 2008.