Исследование вируса птичьего гриппа

Этот пример показывает, как вычислить отношения Ka/Ks для восьми генов в геномах вируса H5N1 и H2N3, и выполнить филогенетический анализ гена HA из вируса H5N1 выделенного из цыплят по всей Африке и Азии. Для филогенетического анализа вы восстановите соседнее дерево и создадите 3-D график расстояний последовательности с помощью многомерного масштабирования. Наконец, вы будете отображать географические местоположения, где каждый последовательность HA был найден на региональной карте. Последовательности, используемые в этом примере, были выбраны из тематического исследования птичьего гриппа на веб-сайте вычислительной геномики [1]. Примечание.В заключительном разделе этого примера требуется Mapping Toolbox™.

Введение

Есть три типа вируса гриппа: Тип A, B и C. Все геномы гриппа состоят из восьми сегментов или генов, которые кодируют для полимеразы B2 (PB2), полимеразы B1 (PB1), полимераза A (PA), гемагглютинин (ХА), нуклеобелок (NP), нейраминидаза (NA), матрица (M1) и неструктурные белки (NS1). Примечание: Вирус типа C имеет гемагглютининэстеразу (HE), гомолог HA.

Из трех типов гриппа тип А может быть наиболее разрушительным. Он поражает птиц (его естественный резервуар), людей и других млекопитающих и был основной причиной глобальных эпидемий гриппа. Тип B поражает только людей, вызывающих локальные эпидемии, и тип C не имеет тенденцию вызывать серьезные заболевания.

Инфлюэнзы типа А далее классифицируются на различные подтипы согласно изменений в аминокислотных последовательностях белков HA (H1-16) и NA (N1-9). Оба белка расположены снаружи вируса. HA присоединяет вирус к хозяину камеры затем помогает в процессе слияния вируса с камерой. NA клипирует недавно созданный вирус из камеры, чтобы он мог перейти к здоровой новой камере. Различие в аминокислотной композиции в белке и рекомбинация различных белков HA и NA способствуют способности инфлюэнза типа A скачкообразным видам хозяина (то есть птице к человеку) и широкой области значений тяжести. Многие новые лекарства разрабатываются для нацеливания на HA и NA белки [2,3,4].

В 1997 году H5N1 подтип вируса птичьего гриппа, вируса гриппа типа А, совершил неожиданный переход к людям в Гонконге, вызвав смерть шести человек. Чтобы контролировать быстро распространяющуюся болезнь, вся птица в Гонконге была уничтожена. Последовательный анализ вируса H5N1 показан здесь [2,4].

Вычисление отношения Ka/Ks для каждого гена H5N1

Исследование коэффициентов Ka/Ks для каждого генного сегмента H5N1 вируса даст некоторое представление о том, как каждый изменяется с течением времени. Ka/Ks - отношение несинонимических изменений к синонимам в последовательности. Для более подробного объяснения коэффициентов Ka/Ks смотрите Анализ Синонимических и несинонимических частот замещения. Чтобы вычислить Ka/Ks, нужна копия гена с двух временных точек. Можно использовать вирус H5N1 выделенный из цыплят в Гонконге в 1997 и 2001 годах. Для сравнения можно включить вирус H2N3 выделенный из кряквы в Альберте в 1977 и 1985 годах [1].

Для целей этого примера данные последовательности предоставляются в четырех структурах MATLAB ®, которые были созданы genbankread.

Загрузите H5N1 и H2N3 данные последовательности.

load('birdflu.mat','chicken1997','chicken2001','mallard1977','mallard1985')

Данные в публичных хранилищах часто кураторятся и обновляются. Вы можете извлечь актуальные наборы данных с помощью getgenbank функция. Обратите внимание, что если данные действительно изменились, результаты этого примера могут несколько отличаться, когда вы используете современные наборы данных.

chicken1997 = arrayfun(@(x)getgenbank(x{:}),{chicken1997.Accession});
chicken2001 = arrayfun(@(x)getgenbank(x{:}),{chicken2001.Accession});
mallard1977 = arrayfun(@(x)getgenbank(x{:}),{mallard1977.Accession});
mallard1985 = arrayfun(@(x)getgenbank(x{:}),{mallard1985.Accession});

Можно извлечь только кодирующий фрагмент нуклеотидных последовательностей, используя featureparse функция. The featureparse функция возвращает структуру с полями, содержащими информацию из раздела Features в файле GenBank, включая поле Sequence, которое содержит только последовательность кодирования.

for ii = 1:numel(chicken1997)
    ntSeq97{ii} = featureparse(chicken1997(ii),'feature','cds','sequence',true);
    ntSeq01{ii} = featureparse(chicken2001(ii),'feature','cds','sequence',true);
    ntSeq77{ii} = featureparse(mallard1977(ii),'feature','cds','sequence',true);
    ntSeq85{ii} = featureparse(mallard1985(ii),'feature','cds','sequence',true);
end

ntSeq97{1}
ans = 

  struct with fields:

                    Location: '<1..>2273'
                     Indices: [1 2273]
    UnknownFeatureBoundaries: 1
                        gene: 'PB2'
                 codon_start: '1'
                     product: 'PB2 protein'
                  protein_id: 'AAF02361.1'
                     db_xref: 'GI:6048850'
                 translation: 'RIKELRDLMSQSRTREILTKTTVDHMAIIKKYTSGRQEKNPALRMKWMMAMKYPITADKRIMEMIPERNEQGQTLWSKTNDAGSDRVMVSPLAVTWWNRNGPTTSTVHYPKVYKTYFEKVERLKHGTFGPVHFRNQVKIRRRVDMNPGHADLSAKEAQDVIMEVVFPNEVGARILTSESQLTITKEKREELKNCNISPLMVAYMLERELVRKTRFLPVAGGTSSVYIEVLHLTQGTCWEQMYTPGGEVRNDDVDQSLIIAARNIVRRATVSADPLASLLEMCHSTQIGGVRMVDILKQNPTEEQAVDICKAAMGLRISSSFSFGGFTFKRTKGFSVKREEEVLTGNLQTLKIQVHEGYEEFTMVGRRATAILRKATRRMIQLIVSGRDEQSIAEAIIVAMVFSQEDCMIKAVRGDLNFVNRANQRLNPMHQLLRHFQKDAKVLFQNWGIEPIDNVMGMIGILPDMTPSTEMSLRGVRVSKMGVDEYSSTERVVVSIDRFLRVRDQQGNVLLSPEEVSETQGMEKLTITYSSSMMWEINGPESVLVNTYQWIIRNWETVKIQWSQEPTMLYNKMEFEPFQSLVPKAARSQYSGFVRTLFQQMRDVLGTFDTVQIIKLLPFAAAPPEQSRMQFSSLTVNVRGSGMRILVRGNSPAFNYNKTTKRLTILGKDAGAITEDPDEGAAGVESAVLRGFLILGKEDKRYGPALSINELSNLTKGEKANVLIGQGDVVLVMKRKRDSSILTDSQTATKRIRMAIN'
                    Sequence: 'agaataaaagaactaagagatttgatgtcgcaatctcgcacacgcgagatactgacaaaaaccactgtggatcatatggccataattaagaagtacacatcaggaagacaggagaagaaccccgctcttagaatgaaatggatgatggcgatgaaatacccgatcacagctgacaaaagaataatggagatgatccctgaaaggaatgagcaaggtcaaactctttggagcaaaacaaatgacgctggatcagacagggtaatggtatcacctctggctgtaacgtggtggaacagaaatggaccaacaacaagtacagtccattatccaaaggtgtataaaacctactttgaaaaggttgaaagattaaaacacggaacctttggccctgttcatttccggaatcaagtcaaaatacgccgcagggttgacatgaaccctggccatgcagatctcagcgctaaagaagcacaagatgtcatcatggaggtcgttttcccaaatgaagttggagccaggatattgacatcagagtcacagctgacaataacaaaggaaaagagggaggaactcaagaattgtaatatttctcctttaatggtggcatatatgttggaaagagaattggttcgcaagaccagattcctaccagtggctggtgggacaagcagcgtatatatagaagtattgcatttgacccaaggaacctgctgggagcagatgtacacaccaggaggggaggtaagaaatgatgatgttgaccaaagtttaatcattgctgctaggaacattgtcaggagagcaacagtatcagcagacccattggcttcactcttggagatgtgccatagcacacaaattggcggagtaagaatggtagacatccttaaacaaaacccaacagaagagcaagctgtagatatatgcaaggcagcaatgggtttaagaatcagctcatccttcagctttggagggttcactttcaaaagaacaaaggggttttctgtcaaaagagaggaagaagtgcttacaggcaacctccaaacattgaagatacaagtacatgaaggatatgaggaattcacaatggttggacgaagagcaacagccattctaagaaaagcaaccagaaggatgatccaactgatagtcagcgggagggacgagcaatcaattgctgaggcaattattgtagcaatggtgttctcacaagaagattgcatgataaaggcagtccgaggtgatttgaatttcgtaaacagagcaaatcaacgactgaaccccatgcaccaactcctgagacacttccaaaaggatgcaaaggtgctgtttcaaaactggggaattgaaccaatcgacaatgtcatggggatgattggaatattgcctgacatgacccccagcacggaaatgtcactaagaggagtgagagttagtaaaatgggggtggatgaatattctagcactgaaagagtggtcgtgagcattgaccgtttcttaagggtccgagatcagcaaggaaatgtactcctatcccctgaagaagttagtgagacacagggaatggaaaagttgacgataacttattcatcgtctatgatgtgggaaattaacggcccagaatcagttctagttaacacataccaatggatcattaggaattgggagactgtaaagatccaatggtcccaagaacccaccatgctatacaataagatggagtttgaaccatttcaatctttagtaccaaaggctgccagaagccaatatagtggatttgtgagaacgctattccagcagatgcgtgatgttttgggaacatttgacactgttcaaataatcaaactactaccatttgcagcagccccacctgaacagagtaggatgcaattttcttctctgactgtgaatgtgaggggatcaggaatgagaatacttgtgagaggtaactcccctgcgtttaactacaacaagacaactaagaggcttacaatacttgggaaggacgcaggtgcaattacagaggacccagatgaaggagcagcaggggtagagtctgcagtattgagagggtttctaattctcggcaaagaagacaaaagatatggaccagcattgagcatcaatgaactgagcaatcttacgaaaggggagaaagctaatgtattgatagggcaaggagacgtggtgttggtaatgaaacggaaacgggactctagcatacttactgacagccagacagcgaccaaaagaattcggatggccatcaatta'

Визуальный осмотр структур последовательности выявил, что некоторые из генов имеют варианты сращивания, представленные в файлах GenBank. Поскольку этот анализ проводится только на PB2, PB1, PA, HA, NP, NA, M1 и NS1 генах, необходимо удалить любые варианты сращивания.

Удалите варианты сращивания из H5N1 1997 года

ntSeq97{7}(1) = [];% M2
ntSeq97{8}(1) = [];% NS2

Удалите варианты сращивания из H2N3 1977 года

ntSeq77{2}(2) = [];% PB1-F2
ntSeq77{7}(1) = [];% M2
ntSeq77{8}(1) = [];% NS2

Удалите варианты сращивания из H2N3 1985 года

ntSeq85{2}(2) = [];% PB1-F2
ntSeq85{7}(1) = [];% M2
ntSeq85{8}(1) = [];% NS2

Вам нужно выровнять нуклеотидные последовательности, чтобы вычислить отношение Ka/Ks. Выравнивание белковых последовательностей для каждого гена (доступно в поле 'translation') с помощью nwalign функцию, затем вставить зазоры в нуклеотидную последовательность с помощью seqinsertgaps. Используйте функцию dnds вычислить несинонимические и синонимичные скорости замещения для каждого из восьми генов в геномах вируса. Если вы заинтересованы в просмотре выравниваний последовательности, задайте для опции 'verbose' значение true при использовании dnds.

Имена генов гриппа

proteins = {'PB2','PB1','PA','HA','NP','NA','M1','NS1'};

H5N1 вирус

for ii = 1:numel(ntSeq97)
    [sc,align] = nwalign(ntSeq97{ii}.translation,ntSeq01{ii}.translation,'alpha','aa');
    ch97seq = seqinsertgaps(ntSeq97{ii}.Sequence,align(1,:));
    ch01seq = seqinsertgaps(ntSeq01{ii}.Sequence,align(3,:));
    [dn,ds] = dnds(ch97seq,ch01seq);
    H5N1.(proteins{ii}) = dn/ds;
end

H2N3 вирус

for ii = 1:numel(ntSeq77)
    [sc,align] = nwalign(ntSeq77{ii}.translation,ntSeq85{ii}.translation,'alpha','aa');
    ch77seq = seqinsertgaps(ntSeq77{ii}.Sequence,align(1,:));
    ch85seq = seqinsertgaps(ntSeq85{ii}.Sequence,align(3,:));
    [dn,ds] = dnds(ch77seq,ch85seq);
    H2N3.(proteins{ii}) = dn/ds;
end
H5N1
H2N3
H5N1 = 

  struct with fields:

    PB2: 0.0226
    PB1: 0.0240
     PA: 0.0307
     HA: 0.0943
     NP: 0.0517
     NA: 0.1015
     M1: 0.0460
    NS1: 0.3010


H2N3 = 

  struct with fields:

    PB2: 0.0048
    PB1: 0.0021
     PA: 0.0089
     HA: 0.0395
     NP: 0.0071
     NA: 0.0559
     M1: 0
    NS1: 0.1954

Примечание: Результаты отношения Ka/Ks могут отличаться от результатов, показанных на [1], из-за вариантов сращивания последовательности.

Визуализируйте коэффициенты Ka/Ks в 3-D гистограммы.

H5N1rates = cellfun(@(x)(H5N1.(x)),proteins);
H2N3rates = cellfun(@(x)(H2N3.(x)),proteins);
bar3([H2N3rates' H5N1rates']);
ax = gca;
ax.XTickLabel = {'H2N3','H5N1'};
ax.YTickLabel = proteins;
zlabel('Ka/Ks');
view(-115,16);
title('Ka/Ks Ratios for H5N1 (Chicken) and H2N3 (Mallard) Viruses');

NS1 HA и NA имеют большие несинонимические и синонимические отношения по сравнению с другими генами как в H5N1, так и в H2N3. Изменения последовательности белков в этих генах были связаны с увеличением патогенности H5N1. В частности, изменения в гене HA могут обеспечить вирусу способность переносить в другие виды, помимо птиц [2,3].

Выполнение филогенетического анализа HA-белка

Вирус H5N1 присоединяется к камерам желудочно-кишечного тракта птиц и дыхательных путей человека. Изменения в белке HA, который помогает связать вирус со здоровой камерой и облегчает его включение в камеру, позволяют вирусу поражать различные органы у одного и того же и разных видов. Это может обеспечить ему способность прыгать от птиц к человеку [2,3]. Вы можете выполнить филогенетический анализ белка HA из вируса H5N1 выделенного из цыплят в разное время (годы) в разных областях Азии и Африки, чтобы исследовать их связь друг с другом.

Загрузите данные аминокислотной последовательности HA из 16 областей/раз из предоставленного MAT-файла birdflu.mat или извлечь актуальные данные последовательности из репозитория NCBI с помощью getgenpept функция.

load('birdflu.mat','HA')
HA = arrayfun(@(x)getgenpept(x{:}),{HA.Accession});

Создайте новый массив структур, содержащий поля, соответствующие аминокислотной последовательности (Sequence) и исходной информации (Header). Можно извлечь исходную информацию из HA с помощью featureparse затем выполните синтаксический анализ с regexp.

for ii = 1:numel(HA)
    source = featureparse(HA(ii),'feature','source');
    strain = regexp(source.strain,'A/[Cc]hicken/(\w+\s*\w*).*/(\d+)','tokens');
    proteinHA(ii).Header = sprintf('%s_%s',char(strain{1}(1)),char(strain{1}(2)));
    proteinHA(ii).Sequence = HA(ii).Sequence;
end

proteinHA(1)
ans = 

  struct with fields:

      Header: 'Nigeria_2006'
    Sequence: 'mekivllfaivslvksdqicigyhannsteqvdtimeknvtvthaqdilekthngklcdldgvkplilrdcsvagwllgnpmcdeflnvpewsyivekinpandlcypgnfndyeelkhllsrinhfekiqiipksswsdheassgvssacpyqgrssffrnvvwlikkdnayptikrsynntnqedllvlwgihhpndaaeqtrlyqnpttyisvgtstlnqrlvpkiatrskvngqsgrmeffwtilkpndainfesngnfiapenaykivkkgdstimkseleygncntkcqtpigainssmpfhnihpltigecpkyvksnrlvlatglrnspqgerrrkkrglfgaiagfieggwqgmvdgwygyhhsneqgsgyaadkestqkaidgvtnkvnsiidkmntqfeavgrefnnlerrienlnkkmedgfldvwtynaellvlmenertldfhdsnvknlydkvrlqlrdnakelgngcfefyhrcdnecmesvrngtydypqyseearlkreeisgvklesigtyqilsiystvasslalaimvaglslwmcsngslqcrici'

Выравнивание аминокислотных последовательностей HA с помощью multialign и визуализируйте выравнивание с seqalignviewer.

alignHA = multialign(proteinHA);
seqalignviewer(alignHA);

Вычислите расстояния между последовательностями с помощью seqpdist с методом Юкса-Кантора. Использование seqneighjoin восстановление филогенетического дерева методом соединения соседей. Seqneighjoin возвращает phytree объект.

distHA = seqpdist(alignHA,'method','Jukes-Cantor','alpha','aa');
HA_NJtree = seqneighjoin(distHA,'equivar',alignHA);

Используйте view метод, сопоставленный с phytree объекты для открытия дерева в Phylogenetic Tree Tool.

view(HA_NJtree);

Визуализация расстояний последовательности с многомерным масштабированием (MDS)

Другой способ визуализации связи между последовательностями - это использование многомерного масштабирования (MDS) с расстояниями, рассчитанными для филогенетического дерева. Эта функциональность обеспечивается cmdscale функция в Statistics and Machine Learning Toolbox™.

[Y,eigvals] = cmdscale(distHA);

Можно использовать собственные значения, возвращенные cmdscale чтобы помочь привести ваше решение о том, использовать ли первые две или три размерности на вашем графике.

sigVecs = [1:3;eigvals(1:3)';eigvals(1:3)'/max(abs(eigvals))];
report = ['Dimension   Eigenvalues    Normalized' ...
          sprintf('\n    %d\t      %1.4f         %1.4f',sigVecs)];
display(report);
report =

    'Dimension   Eigenvalues    Normalized
         1	      0.0062         1.0000
         2	      0.0028         0.4462
         3	      0.0014         0.2209'

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

locations = {proteinHA(:).Header};
figure
plot3(Y(:,1),Y(:,2),Y(:,3),'ok');
text(Y(:,1)+0.002,Y(:,2),Y(:,3)+0.001,locations,'interpreter','no');
title('MDS Plot of HA Sequences');
view(-21,12);

Кластеры, по-видимому, соответствуют группировкам в филогенетическом дереве. Найдите последовательности, принадлежащие каждому кластеру, используя subtree метод phytree. Один из subtreeтребуемыми входами являются номер узла (количество листьев + количество ветвей), который будет корневым узлом нового поддерева. Для вашего примера, кластер, содержащий Хэбэй и Гонконг на графике MDS, эквивалентен поддереву, корневым узлом которого является Branch 14, который является Node 30 (16 оставляет + 14 ветвей).

cluster1 = get(subtree(HA_NJtree,30),'LeafNames');
cluster2 = get(subtree(HA_NJtree,21),'LeafNames');
cluster3 = get(subtree(HA_NJtree,19),'LeafNames');

Получите индекс для последовательностей, принадлежащих каждому кластеру.

[cl1,cl1_ind] = intersect(locations,cluster1);
[cl2,cl2_ind] = intersect(locations,cluster2);
[cl3,cl3_ind] = intersect(locations,cluster3);
[cl4,cl4_ind] = setdiff(locations,{cl1{:} cl2{:} cl3{:}});

Измените цвет и символы маркера на графике MDS, чтобы соответствовать каждому кластеру.

h = plot3(Y(cl1_ind,1),Y(cl1_ind,2),Y(cl1_ind,3),'^',...
    Y(cl2_ind,1),Y(cl2_ind,2),Y(cl2_ind,3),'o',...
    Y(cl3_ind,1),Y(cl3_ind,2),Y(cl3_ind,3),'d',...
    Y(cl4_ind,1),Y(cl4_ind,2),Y(cl4_ind,3),'v');
numClusters = 4;
col = autumn(numClusters);
for i = 1:numClusters
    h(i).MarkerFaceColor = col(i,:);
end
set(h(:),'MarkerEdgeColor','k');
text(Y(:,1)+0.002,Y(:,2),Y(:,3),locations,'interpreter','no');
title('MDS Plot of HA Sequences');
view(-21,12);

Для получения более подробной информации об использовании коэффициентов Ka/Ks, филогенетики и MDS для анализа последовательности см. Cristianini и Hahn [5].

Отображение географических областей вируса H5N1 на карте Африки и Азии

ПРИМЕЧАНИЕ. Вам нужен Mapping Toolbox, чтобы создать следующий рисунок.

Используя инструменты из Mapping Toolbox, можно построить график местоположения, где был выделен каждый вирус, на карте Африки и Азии. Для этого нужны широта и долгота для каждого места. Для получения информации о поиске геопространственных данных в Интернете смотрите Find Geospatial Data Online. Для этого примера предусмотрены широта и долгота для столицы каждой географической области, где были выделены вирусы.

Создайте геоstruct, regionHA, которая содержит географическую информацию для каждого функции или последовательности. Геоstruct должен иметь поля Геометрия, Лат и Лон, которые задают тип функции, широту и долготу. Эта информация используется функциями отображения в Mapping Toolbox для отображения геопространственных данных.

[regionHA(1:16).Geometry] = deal('Point');
[regionHA(:).Lat] = deal(9.10, 34.31, 15.31, 39.00, 39.00, 39.00, 55.26,...
                        15.56, 34.00, 33.14, 34.20, 23.00, 37.35, 44.00,...
                        22.11, 22.11);
[regionHA(:).Lon] = deal(7.10, 69.08, 32.35, 116.00, 116.00, 116.00,...
                        65.18, 105.48, 114.00, 131.36, 131.40, 113.00,...
                        127.00, 127.00, 114.14, 114.14);

Геоstruct может также иметь поля атрибутов, которые содержат дополнительную информацию о каждой функции. Добавьте поля атрибута Name и Cluster к regionHA структура. Поле Cluster содержит номер кластера последовательности, который будет использоваться для идентификации принадлежности кластеров к последовательностям.

[regionHA(:).Name] = deal(proteinHA.Header);

[regionHA(cl1_ind).Cluster] = deal(1);
[regionHA(cl2_ind).Cluster] = deal(2);
[regionHA(cl3_ind).Cluster] = deal(3);
[regionHA(cl4_ind).Cluster] = deal(4);

regionHA(1)
ans = 

  struct with fields:

    Geometry: 'Point'
         Lat: 9.1000
         Lon: 7.1000
        Name: 'Nigeria_2006'
     Cluster: 3

Создайте структуру с помощью makesymbolspec функция, которая будет содержать маркер и спецификации для каждого маркера, отображаемого на карте. Вы передадите эту структуру в geoshow функция. Маркеры и цвета символов соответствуют кластерам на графике MDS.

clusterSymbols = makesymbolspec('Point',...
   {'Cluster',1,'Marker', '^'},...
   {'Cluster',2,'Marker', 'o'},...
   {'Cluster',3,'Marker', 'd'},...
   {'Cluster',4,'Marker', 'v'},...
   {'Cluster',[1 4],'MarkerFaceColor',autumn(4)},...
   {'Default','MarkerSize', 6},...
   {'Default','MarkerEdgeColor','k'});

Загрузите информацию о отображении и используйте geoshow функция для построения графиков расположения вирусов на карте.

load coast
load topo
figure
fig = gcf;
fig.Renderer = 'zbuffer';
worldmap([-45 85],[0 160])
setm(gca,'mapprojection','robinson',...
    'plabellocation',30,'mlabelparallel',-45,'mlabellocation',30)
plotm(lat, long)
geoshow(topo, topolegend, 'DisplayType', 'texturemap')
demcmap(topo)
brighten(.60)

geoshow(regionHA,'SymbolSpec',clusterSymbols);
title('Geographic Locations of HA Sequence in Africa and Asia')

Просмотр географических регионов интереса к Google™ Земли

ПРИМЕЧАНИЕ. Вам нужен Mapping Toolbox, чтобы экспортировать данные в файл в формате KML.

Использование kmlwrite функция из Mapping Toolbox, можно записать информацию о местоположении и аннотации для каждой последовательности в файл в формате KML. Google Earth отображает географические данные из файлов KML в своем браузере Earth. Mapping Toolbox функция переводит геоstruct, такой как regionHA, в файл в формате KML, который будет использоваться Google Earth. Для получения дополнительной информации о kmlwrite см. «Экспорт данных векторной точки в KML».

Можно далее аннотировать каждую последовательность с информацией из раздела Features файла GenBank с помощью featureparse функция. Затем можно использовать эту информацию для заполнения геоstruct, regionHA, и отобразить его в таблицу форме как тег-описание для каждого плацемента в браузере Google Earth. В геоstruct обязательными полями являются поля Геометрия, Лат и Лон. Все другие поля считаются атрибутами знака заполнения.

for i = 1:numel(HA)
    feats = featureparse(HA(i),'Feature','source');
    regionHA(i).Strain = feats.strain;
    if isfield(feats,'country')
        regionHA(i).Country = feats.country;
    else
        regionHA(i).Country = 'N/A';
    end
    year = regexp(regionHA(i).Name,'\d+','match');
    regionHA(i).Year = year{1};
    % Create a link to GenPept record through the accession number
    regionHA(i).AccessionNumber = ...
        ['<a href="http://www.ncbi.nlm.nih.gov/sites/entrez?db=Protein&cmd=search&term=',...
        HA(i).Accession,'">',HA(i).Accession,'</a>'];
end

[regionHA.SequenceLength] = deal(HA.LocusSequenceLength);

Создайте структуру атрибута с помощью makeattribspec функция, которую вы будете использовать для форматирования таблицы описаний для каждого маркера. Структура атрибутов диктует порядок и форматирование каждого атрибута. Можно также использовать его, чтобы не отображать один из атрибутов в геоstruct, regionHA.

attribStruct = makeattribspec(regionHA);

Удалите поле Имя и переупорядочьте поля в структуре атрибутов.

attribStruct = rmfield(attribStruct,'Name');

attribStruct = orderfields(attribStruct,{'AccessionNumber','Strain',...
    'SequenceLength','Country','Year','Cluster'});

regionHA = orderfields(regionHA,{'AccessionNumber','Strain',...
    'SequenceLength','Country','Year','Cluster','Geometry','Lon','Lat',...
    'Name'});

Переформатируйте метки атрибутов для отображения в таблице.

attribStruct.AccessionNumber.AttributeLabel = '<b>Accession Number</b>';
attribStruct.Strain.AttributeLabel = '<b>Viral Strain</b>';
attribStruct.SequenceLength.AttributeLabel = '<b>Sequence Length</b>';
attribStruct.Country.AttributeLabel = '<b>Country of Origin</b>';
attribStruct.Year.AttributeLabel = '<b>Year Isolated</b>';
attribStruct.Cluster.AttributeLabel = '<b>Cluster Membership</b>';

Просмотр файла в Google Earth.

Напишите regionHA геоstruct в файл в формате KML во временной директории.

kmlDirectory = tempdir;
filename = fullfile(kmlDirectory,'HA_geographic_locations.kml');
kmlwrite(filename,regionHA,'Description',attribStruct,'Name',{regionHA.Strain},...
    'Icon','http://maps.google.com/mapfiles/kml/shapes/arrow.png','iconscale',1.5);

Можно отобразить файл KML в браузере Google Earth [6]. В системе должен быть установлен Google Earth. На платформах Windows ® отобразите файл KML с:

winopen(filename)

Для пользователей Unix и MAC отобразите файл KML с:

cmd = 'googleearth ';
fullfilename = fullfile(pwd, filename);
system([cmd fullfilename])

В данном примере файл KML ранее отображался с помощью Google Earth Pro. Изображение Google Earth затем сохранялось с помощью меню Google Earth «Файл - > Сохранить изображение». Так выглядят данные в вашем файле KML при загрузке в Google Earth. Чтобы получить это представление, перемещайтесь и масштабируйте область над Азией.

Щелкните значок, чтобы просмотреть информацию о последовательности. Номер доступа в каждой таблице данных является гиперссылкой на файл последовательности GenPept в базе данных белка NCBI.

Вы можете удалить новый файл KML из выходной директории KML.

delete(filename)
close all

Ссылки

[1] https://computationalgenomics.blogs.bristol.ac.uk/case_studies/birdflu_demo

[2] Laver, W.G., Bischofberger, N. and Webster, R.G., «Disarming Flu Viruses», Scientific American, 280 (1): 78-87, 1999.

[3] Suzuki, Y. and Masatoshi, N. «, Источник и эволюция генов гемагглютинина вируса гриппа», Молекулярная биология и эволюция, 19 (4): 501-9, 2002.

[4] Gambaryan, A., et al. «, Эволюция фенотипа связывания рецепторов вирусов гриппа A (H5)», Вирусология, 344 (2): 432-8, 2006.

[5] Cristianini, N. and Hahn, M.W., «Introduction to Computational Genomics: A Case Studies Approach», Cambridge University Press, 2007.

[6] Изображения Google Earth были получены с помощью Google Earth Pro. Для получения дополнительной информации о Google Earth и Google Earth Pro посетите http://earth.google.com/

Для просмотра документации необходимо авторизоваться на сайте