BLAST

Определение таксономии и функции нуклеотидной последовательности

Ход работы

В качестве нуклеотидной последовательности бралась последовательность из практикума 6
Для поиска использовался blastn с параметрами по умолчанию (megablast), потому что изначально не известно является ли наша последовательность кодирующей
При просмотре GenBank записи первой находки оказалось, что это partial CDS гистона-3
При этом с помощью GenBank было поттверждено, что нуклеотидная последовательность действительно именно кодирующая последовательность гистона-3 (рис.1)

Рис 1. Выдача blastn megablast
Помимо этого был сделан blastx с параметрами по умолчанию (megablast), который ищет находки среди белков, используя нуклеотидную последовательность, это было сделано, чтобы проверить единственность функции этого белка ( т.к. E-value очень маленькое, а Identity, больше 95% у всех лучших находок blastx, можно быть уверенным, что этот ген кодирует высококонсервативную аминокислотную последовательность. Это значит, что очень вероятно моя последовательность выполняет схожую с находками функцию. (рис.2)

Рис 2. Выдача blastx megablast
Для определения близкородственной филогении лучше всего использовать нуклеотидные последовательности, поэтому я взял из лучших находок blastn (рис. 1) разных организмов и построил дл них выравнивание (рис. 3). Хотя подозрительно, что первые 7 лучших находок это гамарусы, но могло так получиться, что из-за меньшего покрытия, чем ,например, у Vibilia, идентичность выше (можно, конечно, попробовать построить дерево и проверить его бутстрепом, но мы этого на тот момент еще не знали. Итого решил взять организмы с E-value < e^-105, потому что с этой степени их значения уменьшались постепенно, а здесь за 13 находок скачок аж на 20 степеней (на следующие 13 находок степень уменьшилась лишь на 4).

Выбранный путь — Рис 3. Выравнивание первых 13 находок blastn (megablast) и моей последовательности 3004

По голому выравниванию сложно что-то сказать вразумительное, видно что в целом последовательность консервативна, а если посмотреть на выравнивание их оттранслированных последовательностей, то видно, что если замены и есть в нуклеотидной последовательности, то большинство из них все же синонимичные (рис. 4)

Поэтому решил построить дерево (рис. 5)

На дереве видно, что последовательность хорошо легла рядом с родом Gammarus (за вид наверное ручаться не могу, хотя и можно предположить, что это Gammarus balcanicus)

Сравнение выдачи разных вариантов запуска blast

Сначала провел разные бласты для моей последовательности megablast, blastn, sensetive blastn
blastn: максимум 1000 находок убрал большой род Niphargus и искал в семействе Gammarida. На выдаче получилось 66 последовательностей (рис. 6)

megablast: максимум 1000 находок убрал большой род Niphargus и искал в семействе Gammarida. На выдаче 25 последовательностей. Что логично, т.к. у blastn длина минимального совпадения 11, а у megablast 28. Yfghbvth megablast не нашел Pseudoniphargus daviui voucher PN0004 histone H3 (H3) gene, partial cds, потому что в нем нету совпадений длинной 28. Так же можно заметить, что порядок последовательностей меняется. Это связано с тем, что уменьшение штрафа за мисматч дает больший вес выравнивания.

sensetive blastn: максимум 1000 находок убрал большой род Niphargus и искал в семействе Gammarida при этом длину слова сделал минимальной 7 и мисматч 1, -1. По хорошему должны ждать находок больше, чем у blastn. Как видно на рис. 9 предсказание сбылось. Однако, теперь к нам попал и мусор с большим E-value, как раз из-за маленькой длины слова и более мягкого мисматча. На рис. 10 продемонстрировано одно из таких выравниваний.

Теперь перейдем к вирусу. Для него выбрал protein=HicA-like toxin gene.
blastn: аналогично, только теперь выкинул всех cellular organisms. На выдаче 132 последовательности.
megablast: как для вирусного blastn. На выдаче 122 последовательности.
sensetive blastn: аналогично для моей последовательности, только выкинул всех cellular organisms. На выходе 139 последовательности, но в них затесались 10 с очень большими E-value. Поэтому можно считать, что их 129.
Общее наблюдение, что выдавались находки исключительно для вирусов Pahexavirus, возможно, это какой-то специфичный для них ген. Судя из названия можно предположить, что это какой-то их токсин.

Поиск белков в геноме

Сначала посмотрим на Amoeboaphelidium protococcarum. В NCBI есть про него пара GenBank записей, из которых берем transl_table=6 для дальнейшего выравнивания.
Искал я белки: ДНК-полимераза дельта - основная полимераза эукариот, гистон H4, и Na+/K+-АТФазу.
Запросы NCBI protein: (gallus gallus[Organism]) AND dna polymerase (здесь взял subunit 3) ; (gallus gallus[Organism]) AND histone h4[Protein Name] ; (Gallus gallus[Organism]) AND ATPase[Protein Name]

В итоге 3 файла и сборка. Далее запускаем сеанс tblastn на кодомо.

makeblastdb -in X5.fasta -dbtype nucl
tblastn -query gallus_dna_pol_d.fasta -db_gencode 6 -db X5.fasta -out dna_pol_d.out
tblastn -query gallus_histone_h4.fasta -db_gencode 6 -db X5.fasta -out histone_h4.out
tblastn -query gallus_atpase.fasta -db_gencode 6 -db X5.fasta -out atpase.out

На выходе имеем для Полимеразы 2 находки. Из них лучшая с весом 31,6 и E-value 1,1 и идентичностью 33%.

> scaffold-17
Length=2125590

 Score = 31.6 bits (70),  Expect = 1.1, Method: Compositional matrix adjust.
 Identities = 14/44 (32%), Positives = 27/44 (61%), Gaps = 2/44 (5%)
 Frame = +3

Query  404      CTDSEDD--FAKTKPPAVPKQPALPVKKEPKEERKNQKKGAATA  445
                CT+S+D   + KT+P A+ +QP+ P  + P ++  N  +  A++
Sbjct  1537893  CTESQDQCTWCKTQPSALEQQPSTPQHRHPLDQIHNSHQICASS  1538024


 Score = 30.8 bits (68),  Expect = 1.8, Method: Compositional matrix adjust.
 Identities = 22/66 (33%), Positives = 35/66 (53%), Gaps = 12/66 (18%)
 Frame = -3

Query  51      ENSGAQLHVTYLV-----------AGNLIQNGHTCHKVAVVREDKLEAMKSKLATVTS-V  98
               +NSG  LH+ YL+           +  L  + HT ++++V     L+  +S LATV + +
Sbjct  231532  QNSGQTLHIEYLIQRVKTIDSNDFSKRLNNSLHTLNQISVSIWFGLQEQQSTLATVLNLI  231353

Query  99      HVYSIQ  104
                VYSIQ
Sbjct  231352  IVYSIQ  231335


> scaffold-105
Length=655906

 Score = 30.4 bits (67),  Expect = 2.2, Method: Compositional matrix adjust.
 Identities = 13/44 (30%), Positives = 28/44 (64%), Gaps = 2/44 (5%)
 Frame = -3

Query  404     CTDSEDD--FAKTKPPAVPKQPALPVKKEPKEERKNQKKGAATA  445
               CT+++D   + KT+P A+ +QP+ P  + P ++  N ++  A++
Sbjct  469775  CTENQDQCTWCKTQPTALEQQPSTPRHRHPLDQIHNSRQICASS  469644

Думаю, что этому не нужно доверять, потому что на мой взгляд уж что-что а субъединица основной полимеразы должна быть очень консервативна. Хотя и второе выравнивание возможно отвечает за 451..458 region_name="PIP-box".
Для АТФазы 33 находки, что возможно было бы не странно, потому что этот белок большой 1000 а.к. Лучшая находка с весом 580, E-value 1e^-177, идентичность 35%

 Score =  580 bits (1495),  Expect = 1e-177, Method: Compositional matrix adjust.
 Identities = 403/1159 (35%), Positives = 594/1159 (51%), Gaps = 171/1159 (15%)
 Frame = +2

Query  8       DKYEPTATSEHGTKKKKAKERDMDELKKEISMDDHKLSLDELHRKY-GTDL-------SR  59
               D+ E +  +    K+   K  D +  KK + + +H LS +EL +KY G          S+
Sbjct  103883  DRVEDSKINTSAVKEPVVKGHDKE--KKSVHITEHMLSPEELSQKYYGVKFDAEKPTQSQ  104056

Query  60      GLTTARAAEILARDGpntltpppttpEWVKFCRQLFGGFSLLLWIGSLLCFLAYGITSVM  119
               GL+ + A E LA  GPN ++PP +TP ++KF   LF  F+ LL    +L F+ YGI    
Sbjct  104057  GLSESLAEEKLAEHGPNMMSPPKSTPAYLKFLHCLFNQFNQLLIFAGILSFILYGI----  104224

Query  120     EGEPNSDNLYLGVVLAAVVIITGCFSYYQEAKSSKIMESFKNMVPQQALVVRNGEKMSIN  179
               + +    +LY+G +L  V +I     + Q  KS  I++SF  +VP+    +R+G+  S+ 
Sbjct  104225  DPQNGISSLYVGAILILVALINSTIEFVQIQKSEAILKSFLGLVPRNCTTIRDGKLKSMG  104404

Query  180     AEGVVVGDLVEVKGGDRIPADLRIISAHGCKVDNSSLTGESEPQTRSPDFSNENPLETRN  239
               AE +V GD++ V+ GD++PAD+ I  A   KVDNSSLTGE++PQ R P  + ENPLE  N
Sbjct  104405  AESLVPGDVIHVRMGDKLPADIYIFWAAEFKVDNSSLTGEADPQERGPGNTQENPLEAHN  104584

Query  240     IAFFSTNCVEGTAVGIVISTGDRTVMGRIASLASGLEGGKTPIAMEIEHFIHLITGVAVF  299
               +AF  +  V G A G+VI TGD TV+G+IA++ASG +  ++P+ +EI++F+ +I  VA  
Sbjct  104585  LAFSGSLAVNGEAYGVVIRTGDFTVLGQIANMASGEKKRQSPMTVEIDNFVKMIALVAAL  104764

Query  300     LGVSFFILSLIL---EYTWLEAVI-----FLIGIIVANVPEGllatvtvcltltaKRMAR  351
                GV FFI+ +     ++   + +I     F IG+ ++ VPE L A VT+ L+  A+RM++
Sbjct  104765  TGVVFFIIGITTLGGQFPDSKQLIAFNFTFAIGVFISWVPEALPAIVTLLLSFAAQRMSK  104944

Query  352     KNCLVKNLEAVGTLGSTSTICSDKTGTLTQNRMTVAHMWFDNQIHEADTTENQSGASFD-  410
               +N LVK+L  V TLGS + + +DKTGTLT+N+MTV + W   ++++   +E++S A+ D 
Sbjct  104945  RNVLVKDLRGVETLGSITLLATDKTGTLTRNQMTVTNFWSSGEMYQVGKSEHKSSANLDG  105124

Query  411     -------------------KSSATWLALSR---------------------------IAG  424
                                  + +   L L+R                           I  
Sbjct  105125  QPKSAGDENMPPVINGKDGEQAIGMLPLNRVPTVQMAPADDFKQFSLDVKSLNDLSMICY  105304

Query  425     LCNRAVFQANQENVPILKRAVAGDASESALLKCIELCCGSVKEMRERYPKVVEIPFNSTN  484
               LC++A F A   +VP+ +R+V GDA+ES L             + E+YPKV EIPFNSTN
Sbjct  105305  LCSKARFDAT--DVPLKQRSVIGDATESGLFLFAANTLPDSDTLAEQYPKVFEIPFNSTN  105478

Query  485     KYQLSIHKNANAGESRHLLVMKGAPERILDRCDSILIH-GKVQPLDEEIKDAFQNAYlel  543
               K+ LSIHK  +      LL+ KGAPERI   C  I    G  + +  E+K  FQ +Y  L
Sbjct  105479  KWHLSIHKMKHDNGDLTLLI-KGAPERIFRLCSKIYSQSGDEKEITAEMKQDFQKSYEAL  105655

Query  544     gglgervlgFCHLALPDDQFPEGFQFDTDEVN--------------FPVEKLCFVGLMSM  589
                  G RV+G     LP   +P  F F  D+ +              +P     F GL+S+
Sbjct  105656  AAKGHRVIGTAKFNLPAADYPADFTFRKDDQDKKELGEDKPGILGTYPKGGYTFCGLVSL  105835

Query  590     IDPPRAAVPDAVGKCRSAGIKVIMVTGDHPITAKAIAKGVGIISDGNETVEDIAARLNIP  649
                DPP+  V +A+GKCR AG+KV+MVTGDHP+TA+AIA+ + ++    ET E +A R N  
Sbjct  105836  EDPPKHGVREAIGKCRQAGVKVMMVTGDHPLTAEAIARKINLML--QETKEMVAKRTNRA  106009

Query  650     VSQVNPRDAKACVVHGSDLKDMTSEQLDDILLHHTEIVFARTSPQQKLIIVEGCQRQGAI  709
               +S +   +  + V+HG  +  +T +  + IL +  EI+FARTSP+ KL IV+ CQ +G I
Sbjct  106010  LSSIQEHEYNSIVIHGEKVDSLTEDDWERIL-NKDEIIFARTSPKHKLQIVKHCQERGHI  106186

Query  710     VAVTGDGVNDSPALKKADIGVAMGIAGSDVSKQAADMILLDDNFASIVTGVEEGRLIFDN  769
               V VTGDGVNDSPALKKAD+G++M I+GSDVSK+AA MILLDDNFAS V G+ EGRLIF N
Sbjct  106187  VGVTGDGVNDSPALKKADLGISMNISGSDVSKEAAAMILLDDNFASTVHGISEGRLIFQN  106366

Query  770     LKKSIAYTLTSNIPEITPFLIFIIANIPLPLGTCTILCIDLGTDMVPAISLAY---EQAE  826
               LKK + YTL   +PE+   L+F++  IPLPL    I+ IDLG +    +S A+   E  E
Sbjct  106367  LKKCVRYTLCHILPEVIANLLFVVVPIPLPLYALQIILIDLGFEFFNGLSYAWEVPEHGE  106546

Query  827     SDIMKRQPRNPKTDKLVN-------------ERLISMAYGQIGMIQALGGFFTYFVIMAE  873
                 I+   PR P + + ++              +++ M YG      A        V  +E
Sbjct  106547  EGILMALPRKPVSLRSIDLLRRNNAYKAKQQNQVMQMVYGNNRPAAAAENSLVQDVDTSE  106726

Query  874     NGFLPSGLVGIRLQWDDRWINDVEDSYGQQWTFEQRK--------IVEFT----------  915
                 +  S  +  +      W   + D  G Q  +E R+        I+ +           
Sbjct  106727  MNW--STYLKYKTAATRVWFRSLFDKEGWQLYWEPRQEETLVDADILSYAYLEVGVLETI  106900

Query  916     -CHTAFF-VSIVVVQWADLIICKT----RRNSVFQQGMKNKIL--------------IFG  955
                C  AFF V  V + W   ++ +        +   Q +K  +               IF 
Sbjct  106901  GCLLAFFHVFYVELGWTPAVVAQNAAQFSNTNALPQNLKEALFKAQSAYYFALLIMQIFN  107080

Query  956     LF-------------------------EETALAAFLSYCPGMDVALRMYPLKPTWWFCAF  990
               +F                            A + F+ Y P M    +   L P +W   F
Sbjct  107081  MFCCKVTTSYPFGWRVLKNKVTWISLAVSVAFSCFVIYPPFMHDVFQTNYLSPQYWLFPF  107260

Query  991     PYSLLIFLYDEIRKLIIRR  1009
                   +I +Y   R    RR
Sbjct  107261  VMGFVILMYVSARVAYRRR  107317

Опять же я бы не поверил, что АТФ-аза может так различаться, хот и E-value многообещающий и при этом последовательность целиком легла на скэффолд, что тоже подозрительно.
Для гистона H4 все получше потому что лучшая выборка имеет массу 159, E-value 6e^-49 и при этом 94% идентичность.

 Score =  159 bits (403),  Expect = 6e-49, Method: Compositional matrix adjust.
 Identities = 77/82 (94%), Positives = 81/82 (99%), Gaps = 0/82 (0%)
 Frame = -2

Query  22   VLRDNIQGITKPAIRRLARRGGVKRISGLIYEETRGVLKVFLENVIRDAVTYTEHAKRKT  81
            +LRDNIQGITKPAIRRLARRGGVKRISGLIYEETRGVLK FLENV+RD+VTYTEHAKRKT
Sbjct  331  ILRDNIQGITKPAIRRLARRGGVKRISGLIYEETRGVLKSFLENVVRDSVTYTEHAKRKT  152

Query  82   VTAMDVVYALKRQGRTLYGFGG  103
            VTA+DVVYALKRQGRTLYGFGG
Sbjct  151  VTALDVVYALKRQGRTLYGFGG  86

Вот здесь впринципе можно быть уверенным, что нашлось совпадение.

Поиск гена в контиге

Контиг взял из сборки с AC: GCA_000767325.1 для Danio rerio. Длина 22667 н.о.
Далее провожу blastx, который ищет находки среди белков по нуклеотидной последовательности. Для меньшего количества ненужных находок будем искать только у эукариот и только по проверенным reference белкам (соответствующая бд выбрана: refseq_protein)
Итого на выдаче 64 последовательности рис. 11

Теперь отбросим все находки с E-value >= 0,01 и среди них обратим внимание на находки с Ident, хотя бы больше 70%, то вероятно можем предположить, что эти гены этих белков расположены на выбранном контиге