Идентификация потенциальных ключевых генов, являющихся маркерами стресса засухи в листьях сои
Идентификация потенциальных ключевых генов, являющихся маркерами стресса засухи в листьях сои
Аннотация
Соя, одна из важнейших и наиболее распространенных бобовых культур в мире, подвергается абиотическим стрессам, таким как засуха и засоление почвы, которые снижают рост и урожайность. Анализ данных РНК-секвенирования листьев сои из открытых баз данных проводится с целью выявить гены-маркеры стресса в условиях засухи и засоления почвы. С помощью биоинформатического анализа сформирован список генов-маркеров стресса, наиболее сильно изменивших регуляцию. Среди наиболее повысивших экспрессию генов есть такие, которые участвуют в регуляции роста, развития растения и реакции на абиотический стресс. Наиболее сильно подавилась экспрессия генов, связанных с процессами, требующими затрат энергии, а также участвующих в росте и развитии.
1. Введение
Стрессовые условия окружающей среды оказывают влияние на рост и развитие сельскохозяйственных культур, вызывая травмы и повреждения растений, которые могут привести к снижению урожайности . По сравнению с рекордными урожаями абиотические стрессы, такие как засуха, засоление почвы, холод и жара, могут снизить урожайность в среднем более чем на 60% . Одними из абиотических стрессовых факторов для растений являются засуха и засоление. Из-за изменения климата и давления со стороны растущего населения Земли эти стрессы становятся все более серьёзными . Также засуха и засоление почвы ограничивают площади выращивания большинства культур .
Соя (Glycine max L.), одна из важнейших бобовых культур, является богатым источником белков, ненасыщенных жиров, углеводов и клетчатки, используются в производстве продуктов питания, кормов для животных и сырья для промышленности, а также применяется в биотехнологии , .
Рост и развитие сои очень чувствительны к дефициту воды, а также к засоленной почве , . Качество и урожайность сои сильно страдают от засухи . Среди сельскохозяйственных культур соя обладает наибольшей чувствительностью к засухе. Засухоустойчивость — важный признак, связанный с урожайностью. Засухоустойчивым считается тот сорт сои, урожайность которого намного выше, чем у других сортов в условиях засухи . Влияние засухи на урожайность также зависит от стадии развития растения, на которой оно подвергается воздействию засухи. Растения, подвергшиеся засухе на репродуктивной стадии, более уязвимы к дефициту воды . На этой стадии потребность сои в воде увеличивается вдвое по сравнению с вегетативной . Продолжительность, время и интенсивность засухи играют ключевую роль в определении урожайности и качества сои. При постоянном воздействии засухи растения могут испытывать окислительный стресс, повреждение мембранных структур, потерю ионов из клеток и денатурацию белков .
Солевой стресс также отрицательно влияет на рост, развитие и урожайность растений. Сильно засоленная почва может образоваться в результате как естественных процессов, например, выветривания или дождей с содержанием морской соли, так и деятельности человека — орошения, расчистки земель . Соя считается чувствительной к соли культурой. Токсичное накопление Na+ и Cl− нарушает нормальный процесс фотосинтеза, повреждая рецепторную часть фотосистемы II в листьях сои, приводит к окислительному стрессу, который негативно сказывается на росте и метаболизме растений, причем Cl- для проростков сои более токсичен, чем Na+ , .
В течение жизненного цикла растения подвергаются воздействию различных абиотических стрессовых факторов, влияющих на их рост и развитие. Чтобы адаптироваться к изменениям окружающей среды, растения выработали определенные механизмы регуляции. В ответ на стрессовые сигналы индуцируется экспрессия генов, связанных с устойчивостью к стрессу, что приводит к метаболическим и физиологическим изменениям, которые обеспечивают адаптацию растений к новым условиям окружающей среды , . У разных генотипов одного и того же вида растений способность переносить абиотические стрессы значительно различается. Для выявления генов, участвующих в формировании устойчивости, проводятся исследования регуляторных сетей, отвечающих за реакцию на стресс .
Понимание механизмов реагирования на абиотический стресс важно не только для изучения физиологии растений, но и для разработки потенциальных способов улучшения сельскохозяйственных культур . Достижения в области высокопроизводительных технологий, таких как секвенирование РНК и микрочипы, сделали доступным большое количество транскриптомных данных. В результате есть возможность получать более надежные результаты, используя информацию из нескольких источников. Кроме того, можно изучать экспрессию тысяч генов при различных абиотических стрессах .
Обработка данных секвенирования транскриптома, сделанная для разных растений, помогает выявить регулируемые гены и предоставляет обширный ресурс для улучшения растений. Такие результаты дают представление о генах, реагирующих на абиотический стресс, и предлагают новые целевые гены для инженерии устойчивости к нему и для растений сои . Число анализов транскриптома резко увеличилось в последнее время, и в большинстве таких исследований проводилось секвенирование РНК (RNA-Seq) с использованием платформы Illumina. Такое обилие исследований делает возможным метаанализ данных RNA-Seq с использованием общедоступных данных .
В этом исследовании проводится анализ данных РНК-секвенирования из открытых баз данных в условиях засухи и засоления, с целью выявить гены-маркеры стресса растений сои в этих условиях.
2. Подготовка сырых данных
Условия засухи и засоления выбраны потому, что в ответ на них активируется модуль MAPK, функционирующий в сигнальной системе абсцизовой кислоты. Абсцизовая кислота опосредует реакции на стресс окружающей среды, главным образом на водный , , . Из-за засоления почвы, то есть низкого качества воды, у растений также возникает стресс засухи .
Данные RNA-Seq взяты из общедоступной базы данных Национального центра биотехнологической информации Gene Expression Omnibus (NCBI GEO) . Комплексный поиск в GEO был выполнен с использованием поискового запроса: для засоления ((«Glycine max»[Organism] OR Glycine max[All Fields]) AND («sodium chloride»[MeSH Terms] OR salt[All Fields])) AND «Glycine max»[porgn], для засухи («Glycine max»[Organism] OR Glycine max[All Fields]) AND («water»[MeSH Terms] OR Water[All Fields]) AND Deficit[All Fields] и («Glycine max»[Organism] OR Glycine max[All Fields]) AND Drought[All Fields]. Тип исследования — профилирование экспрессии с помощью высокопроизводительного секвенирования (Expression profiling by high throughput sequencing). Всего 2 исследования, показанные в таблице 1, во обоих парноконцевые данные.
Во обоих исследованиях отбирались листья. Условия засухи были во всех, соленость — только в одном из них. Количество контрольных образцов и обработанных стрессом составило 6 и 9 соответственно.
Таблица 1 - Исследования, данные SRR из которых использованы для обработки
Исследование | Образец, контроль | Образец, экспериментальный | Воздействие (обработка) |
SRP450502 | SRR25360934 | SRR25360931 | засуха (30% воды, доступной для транспирации), соленость (15 mM NaCl) |
SRR25360932 | |||
SRR25360935 | SRR25360933 | ||
SRR25360922 | |||
SRR25360936 | SRR25360923 | ||
SRR25360924 | |||
SRP058975 | SRR2048194 | SRR2048191 | засуха 24 ч |
SRR2048195 | SRR2048192 | ||
SRR2048196 | SRR2048193 |
В исследовании SRP450502 дефицит воды создавался через 21 день после начала эксперимента, а листья были собраны с 5–7 разных растений через 10 дней после начала дефицита воды. Ткани растений сои, подвергшихся воздействию дефицита воды, были возрастом 31 день. В SRP058975 измерения проводились в течение 6 периодов суток, на стадии развития растений V1, в таблице 1 взяты данные только в один период (24 часа).
Сырые данные загружали с ncbi с помощью инструмента SRA-toolkit (версия 3.2.0) с использованием команды prefetch --option-file. Для извлечения данных из архивов sequence read archive (SRA), преобразования в формат fastq, сжатия в формат gz и обработки парных чтений применялась команда fastq-dump --gzip --split-3 SRR*/*.sra. Качество прочтений оценивалось с использованием FastQC (0.12.0), обобщение результатов — с помощью MultiQC (1.27.1) , . Для поиска и удаления нежелательных последовательностей из считываний (адаптеров и др.) и обработки последовательностей плохого качества использовались программы Trimmomatic (версия 0.39) и BBTools (37.62) , . Для количественной оценки транскриптов использовалась программа Salmon (v.1.10.3) (индексирование и количественная оценка), сопоставлялось с эталонной последовательностью кДНК, загруженной с сайта Ensembl Plants , . Для нормализации и выявления дифференциально экспрессируемых генов использовали pyDESeq2 (v.0.5.0) . Достоверно дифференциально экспрессируемыми генами считались гены со значениями p-value меньше 0,05, их число составило 759.
По умолчанию DESeq2 использует метод Бенджамина-Хохберга для коррекции p-value. С учетом скорректированных этим методом значений p-value количество достоверно экспрессируемых генов снижается с 759 до 3 при условии, если принять p-value не больше 0,05 и до 7 генов, если p-value не более 0,1. Такое сильное снижение количества достоверных генов после применения FDR-коррекции показывает недостаток мощности данных и потерю важных биологических сигналов. Эта коррекция сильно ужесточила критерии значимости, препятствуя полноценному изучению процессов умеренной экспрессии. Вместо этого был проведен первичный анализ без коррекции с обычными p-value, позволяющий сохранить максимальное число генов-кандидатов для последующего этапа проверки и анализа.
С помощью скрипта на основе данных полученной таблицы pyDESeq2 получили график с осями -log10p-value и log2foldchange, на котором показаны достоверно и недостоверно экспрессируемые гены. Затем проведен анализ обогащения генной онтологии (ГО) на сайте ShinyGO для определения основных биологических процессов, с которыми связаны выявленные дифференциально экспрессируемые гены .
3. Основные результаты
3.1. Оценка качества
Суммарный начальный размер архивов до обработки был 13,4 Гб для контрольных условий и 22,5 Гб для экспериментальных. Обрабатывали парноконцевые данные. Во всех образцах из обеих групп проводили обрезку рядов программой Trimmomatic. После обработки Trimmomatic среди выходных файлов далее работали только с парами чтений (_1P и _2P), не использовали одиночные чтения (_1U и _2U).
Некоторые образцы (SRR2048195_2P, SRR2048195_2P и SRR2048192_2P), качество которых осталось недостаточным, триммировали с помощью BBTools.
В таблице 2 показано количество последовательностей (в миллионах (M seqs), это сумма Seq для обоих считываний) в образцах до триммирования и после. Так как данные парноконцевые, то для каждого образца есть два варианта: _1 и _2.
Таблица 2 - Количество исходных и отфильтрованных последовательностей после триммирования для образцов
Образец | Количество последовательностей до обрезки (M seqs) | Количество последовательностей после обрезки (M seqs) |
SRR2048194_1 | 20,1 | 15,9 |
SRR2048194_2 | 20,1 | 15,9 |
SRR2048195_1 | 20,1 | 14,4 |
SRR2048195_2 | 20,1 | 7,4 |
SRR2048196_1 | 47,7 | 35,5 |
SRR2048196_2 | 47,7 | 18,2 |
SRR25360934_1 | 22,8 | 22,2 |
SRR25360934_2 | 22,8 | 22,2 |
SRR25360935_1 | 29,2 | 28,4 |
SRR25360935_2 | 29,2 | 28,4 |
SRR25360936_1 | 26,9 | 26,4 |
SRR25360936_2 | 26,9 | 26,4 |
SRR2048191_1 | 44,3 | 33,5 |
SRR2048191_2 | 44,3 | 33,5 |
SRR2048192_1 | 49,4 | 33,5 |
SRR2048192_2 | 49,4 | 17,7 |
SRR2048193_1 | 52,4 | 42,1 |
SRR2048193_2 | 52,4 | 42,1 |
SRR25360922_1 | 21,4 | 20,7 |
SRR25360922_2 | 21,4 | 20,7 |
SRR25360923_1 | 23,0 | 22,4 |
SRR25360923_2 | 23,0 | 22,4 |
SRR25360924_1 | 20,5 | 19,9 |
SRR25360924_2 | 20,5 | 19,9 |
SRR25360931_1 | 23,6 | 22,9 |
SRR25360931_2 | 23,6 | 22,9 |
SRR25360932_1 | 22,5 | 21,7 |
SRR25360932_2 | 22,5 | 21,7 |
SRR25360933_1 | 24,1 | 23,4 |
SRR25360933_2 | 24,1 | 23,4 |
Из таблицы можно заметить, что после триммирования количество последовательностей уменьшилось, так как были удалены адапторные и другие технические последовательности. Это важно, поскольку эти последовательности могут мешать последующему анализу и снижают качество.
3.2. Дифференциальная экспрессия
На основе данных полученной таблицы pyDESeq2 с помощью скрипта получили график «вулкан» с осями -log10P-value и log2foldchange (рисунок 1), на котором показаны точками достоверно и недостоверно экспрессируемые гены и их количество. Среди достоверно экспрессируемых генов up-reg — гены, экспрессия которых повысилась, down-reg — гены, экспрессия которых понизилась. Достоверными считались точки со значениями p-value больше 0,5.
Рисунок 1 - Volcano plot дифференциально экспрессируемых генов:
up-reg - гены, увеличившие экспрессию; down-reg - гены, снизившие экспрессию
Рассмотрим 25 генов с наиболее повышенной экспрессией, список которых показан в таблице 3, и 25 с наиболее пониженной, показанных в таблице 4. Выбраны на основе параметра log2FoldChange, который показывает насколько сильно изменилась экспрессия генов в экспериментальных условиях по сравнению с контрольными. Для увеличивших экспрессию генов список сортировали по наибольшим значениям log2FoldChange, для уменьшивших — наоборот, по наименьшим значениям этого параметра. Описание генов взято с сайтов NCBI и UniProt . Названия генов Glyma определяли через индексы транскриптов KRH, KRG и RCW с помощью сайта Ensembl Plants.
Таблица 3 - 25 генов, экспрессия которых наиболее сильно увеличилась
Ген | log2FoldChange | Описание | Функция |
GLYMA_13G021800 | 10,1150765 | hypothetical protein | N/A |
GLYMA_13G311000 | 9,475082861 | fasciclin-like arabinogalactan protein 15 | N/A |
GLYMA_11G175700 | 9,374786653 | E3 ubiquitin protein ligase | ubiquitin protein ligase activity/chromatin organization/protein ubiquitination |
GLYMA_09G223300 | 8,931801879 | DNA topoisomerase 1 beta | chromosome segregation/DNA topological change |
GLYMA_08G264900 | 8,899239799 | Phototropin-2 (non-specific serine/threonine protein kinase) | blue light photoreceptor activity |
GLYMA_07G261400 | 8,646048704 | lactoylglutathione lyase | methylglyoxal catabolic process to D-lactate via S-lactoyl-glutathione |
GLYMA_17G102800 | 8,516915567 | Formin-like protein | phosphoprotein phosphatase activity |
GLYMA_05G168100 | 8,41744794 | bZIP transcription factor bZIP100 | DNA-binding transcription factor activity |
GLYMA_17G005100 | 8,336940928 | phosphatidylinositol/phosphatidylcholine transfer protein SFH2 (CRAL-TRIO domain-containing protein) | N/A |
GLYMA_07G223400 | 8,17456152 | O-fucosyltransferase 23 | fucose metabolic process/ pollen-pistil interaction |
GLYMA_04G230800 | 8,173751363 | uncharacterized LOC100792159 | N/A |
GLYMA_18G136900 | 7,991056092 | serine/arginine-rich splicing factor RSZ21-like | mRNA cis splicing, via spliceosome |
GLYMA_08G331200 | 7,989181625 | ER membrane protein complex subunit 3 | N/A |
GLYMA_10G239500 | 7,946552525 | valine--tRNA ligase, mitochondrial 1 | valyl-tRNA aminoacylation / post-embryonic development / reproductive structure development |
GLYMA_03G052800 | 7,863521193 | TIR domain-containing protein | S-adenosylmethionine-dependent methyltransferase activity / defense response / signal transduction |
GLYMA_04G191600 | 7,836373058 | UDP-N-acetylglucosamine--peptide N-acetylglucosaminyltransferase SPINDLY | gibberellic acid mediated signaling pathway |
GLYMA_10G266900 | 7,766813625 | sarcoplasmic reticulum histidine-rich calcium-binding protein | RNA splicing |
GLYMA_15G101900 | 7,696307644 | putative cobalt ion transmembrane transporter (ABC transmembrane type-1 domain-containing protein) | cobalt ion transmembrane transporter activity |
GLYMA_20G009200 | 7,661347024 | hypothetical protein | N/A |
GLYMA_10G245400 | 7,646573831 | NAB domain-containing protein | actin filament binding |
GLYMA_11G086000 | 7,646573831 | DUF4378 domain-containing protein | N/A |
GLYMA_08G013200 | 7,611169647 | probable LRR receptor-like serine/threonine-protein kinase At1g63430 (Protein kinase domain-containing protein) | protein kinase activity |
GLYMA_05G046000 | 7,596840017 | cyclin-dependent kinase | protein serine/threonine kinase activity |
GLYMA_08G054100 | 7,592053918 | ADP/ATP translocase | ATP:ADP antiporter activity (Catalyzes the exchange of ADP and ATP across the membrane) |
GLYMA_07G009900 | 7,492522691 | Glucose-methanol-choline oxidoreductase N-terminal domain-containing protein | oxidoreductase activity, acting on CH-OH group of donors |
В таблице 4, наоборот, показаны 25 генов, экспрессия которых наиболее сильно уменьшилась, сортировка по наименьшим значениям параметра log2FoldChange.
Таблица 4 - 25 генов, экспрессия которых наиболее сильно уменьшилась
Ген | log2FoldChange | Описание | Функция |
GLYMA_11G084700 | -9,666622991 | Uncharacterized protein | полиубиквитин-зависимое связывание белков / ubiquitin-dependent protein catabolic process |
GLYMA_02G136800 | -8,812982149 | PX domain-containing protein | phosphatidylinositol binding |
GLYMA_16G041200 | -8,376361571 | glutamate dehydrogenase 1 | glutamate catabolic process |
GLYMA_02G208700 | -8,28958539 | 14-3-3 protein SGF14g | Взаимодействие с фосфосерином на взаимодействующем белке |
GLYMA_08G086100 | -8,218512686 | Protein CHROMATIN REMODELING 4 | ATP-dependent chromatin remodeler activity |
GLYMA_11G106300 | -8,168778873 | protein NUCLEAR FUSION DEFECTIVE 4 (Nodulin-like domain-containing protein) | transmembrane transporter activity |
GLYMA_19G037200 | -8,168778873 | protein furry homolog (Cell morphogenesis protein N-terminal domain-containing protein) | cell morphogenesis |
GLYMA_15G132500 | -8,159017331 | Fatty acid hydroxylase domain-containing protein | fatty acid biosynthetic process |
GLYMA_20G227200 | -8,069311246 | protein CHROMATIN REMODELING 20 (ATP-dependent helicase ATRX) | ATP-dependent chromatin remodeler activity / transcription coregulator activity |
GLYMA_11G034600 | -8,069311246 | DNA mismatch repair proteins mutS family domain-containing protein | ATP-dependent DNA damage sensor activity / meiotic cell cycle |
GLYMA_12G007500 | -8,069311246 | PAP36 purple acid phosphatase 36 | acid phosphatase activity |
GLYMA_18G199400 | -7,970367052 | RNA polymerase II C-terminal domain phosphatase-like 2 (protein-serine/threonine phosphatase) | histone H2AXS140 phosphatase activity / RNA binding |
GLYMA_16G005600 | -7,958432482 | DNA topoisomerase 2-binding protein 1-A (BRCT domain-containing protein) | N/A |
GLYMA_07G261400 | -7,948586317 | lactoylglutathione lyase | lactoylglutathione lyase activity |
GLYMA_19G071800 | -7,878123882 | Sulfite oxidase | oxidoreductase activity |
GLYMA_08G044700 | -7,816925923 | zinc finger protein BRUTUS | inorganic ion homeostasis |
GLYMA_15G230300 | -7,816895562 | retinoblastoma-related protein 1 | regulation of cell cycle / regulation of transcription by RNA polymerase II |
GLYMA_16G220100 | -7,757947948 | phospholipid-transporting ATPase 1 | ATPase-coupled intramembrane lipid transporter activity |
GLYMA_15G166700 | -7,72827923 | DNA ligase IV | obsolete DNA ligation involved in DNA repair / DNA recombination |
GLYMA_10G238100 | -7,707739632 | homeobox-leucine zipper protein ANTHOCYANINLESS 2 | DNA-binding transcription factor activity, RNA polymerase II-specific |
GLYMA_12G059800 | -7,628601473 | transcription initiation factor IIA large subunit (Transcription factor IIA alpha/beta subunit) | transcription initiation at RNA polymerase II promoter |
GLYMA_08G044700 | -7,578469063 | CHY-type domain-containing protein | inorganic ion homeostasis, ubiquitin-dependent protein catabolic process |
GLYMA_08G081600 | -7,558668068 | HEN1A protein Hua enhancer 1a (Small RNA 2'-O-methyltransferase) | RNA methylation / siRNA processing |
GLYMA_07G127800 | -7,529450793 | IST1-like protein | protein transport |
GLYMA_10G238500 | -7,46955387 | Protein FAR1-RELATED SEQUENCE | regulation of DNA-templated transcription |
3.3. Построение обогащения онтологии генов
Проведено построение обогащения генной онтологии на сайте ShinyGO для того, чтобы определить какие биологические процессы участвовали в реакции на стресс. Использовали список из достоверно экспрессируемых генов. Из 759 идентификаторов 725 (96%) сопоставлены с генами glycine max. Определен 21 биологический процесс GO. Параметр FDR (False Discovery Rate) 0,05. С помощью этого параметра оценивается статистическая значимость. База данных путей — биологический процесс GO. На рисунке 2 показана диаграмма биологический процесс GO/кратное обогащение, построенная в ShinyGO на основе списка генов, изменивших регуляцию.
Рисунок 2 - Диаграмма биологический процесс GO/кратное обогащение
Из-за того, что получили неожиданные термины GO для листьев, также проведен анализ обогащения на сайте AgriGO . Использовано (аннотировано) только 373 гена из задаваемого списка, число значащих GO терминов 21, уровень значимости 0.05, статистический метод тестирования - критерий Фишера. Результаты анализа в виде иерархического графа показаны на рисунке 3.
Рисунок 3 - Графические результаты обогащения GO для каждой из трех категорий:
а – биологический процесс; б – клеточный компонент; в – молекулярная функция
4. Обсуждение
При абиотических стрессах активируется множество различных механизмов. Они варьируются от плазмолиза до окислительного стресса, от нарушения потока электронов до биомолекул, клеточной деградации/дисфункции других. Транскриптом, протеом и метаболом растения модулируется для защиты от стресса и восстановления . В растениях выработалось множество защитных механизмов для регуляции экспрессии определенных генов с помощью различных сигнальных путей, изменения использования воды, энергетического обмена и поглощения АФК, чтобы адаптироваться к неблагоприятным условиям окружающей среды .
4.1. Повышение
Среди генов, наиболее повысивших экспрессию несколько таких, которые участвуют в регуляции роста, развития растений и их реакции на абиотический стресс. Среди них формин-подобный белок, фасциклиноподобный арабиногалактановый белок, которые участвуют в адаптации растений к стрессу , . Формин-подобный белок и белок, содержащий домен NAB участвуют в регуляции актина, контролируя структуру и динамику актинового цитоскелета , . Повысилась экспрессия гена убиквитинлигазы E3, участвующего в убиквитин-опосредованной деградации белков, которая является важным процессом клеточной защиты от воздействия окружающей среды . Есть гены, участвующие в сигнальных путях абсцизовой кислоты, связанных с реакцией на водный стресс. Некоторые из них — транскрипционные факторы, такие как bZIP, которые имеют ключевую роль в реакции на стресс . Они регулируют экспрессию нижестоящих генов, отвечающих за устойчивость растений к стрессу. Повысился ген циклин-зависимой киназы, который играет важную роль в регуляции клеточного цикла и имеет стресс-реактивную природу. В условиях биотического и абиотического стресса на рост растений негативно влияет ингибирование клеточного цикла из-за удлинения фазы S и задержки вступления в митоз. В других работах также были обнаружены его высокие уровни экспрессии в листьях и в корнях при засухе, обработке абсцизовой кислотой и солевом стрессе, по сравнению с необработанными стрессом растениями . Кроме того, среди активированных генов есть несколько генов гипотетических белков, в том числе наиболее сильно активированный ген. Повысилась экспрессия гена белка, содержащего TIR-домен. Домены TIR — это ферменты, катализирующие выработку малых молекул для передачи иммунных сигналов и приводящих к гибели растительных клеток . Активировался ген возможной LRR рецептороподобной серин треониновой протеин киназы. Рецептороподобные киназы играют центральную роль в передаче сигналов во время распознавания патогена, последующей активации защитных механизмов растений и контроля развития . Повысилась экспрессия богатого гистидином кальций-связывающего белка саркоплазматического ретикулума. Он связывает Ca(2+). Изменение концентрации кальция в цитозоле является одной из первых реакций на различные раздражители, и элементы, переносящие кальций, активно поддерживают этот поток и гомеостаз . Повысилась экспрессия фототропина-2, чувствительного к ультрафиолетовому и синему свету фоторецептора растений. Он контролируют перемещение ядер, избегание хлоропластов и тёмное положение .
Увеличилась экспрессия O-фукозилтрансферазы, ответственного за посттрансляционные модификации (O-гликозилирование ядерных и цитозольных белков), которые динамически регулируют различные физиологические и морфологические процессы на протяжении всей жизни растений, от развития эмбриона до развития цветка . Также повысилась экспрессия трансмембранного ABC-транспортера ионов кобальта. ABC (ATP-binding cassette)-транспортерные системы имеют решающее значение для роста и развития растений . Повысилась экспрессия белка с доменом неизвестной функции (DUF). Это белковый домен, функция которого не определена.
4.2. Понижение
Из генов, понизивших экспрессию, есть такие, которые связаны с процессами, требующими затрат энергии, например, фосфолипид-транспортирующая АТФаза 1, также связанные с ростом и развитием, так как во время стресса замедляется рост и развитие растения, поэтому они подавились. Также снизился ген АТФ-зависимой геликазы ATRX. ATRX подвергается зависимому от клеточного цикла фосфорилированию, которое регулирует его связь с ядерным матриксом и хроматином и предполагает его участие в регуляции генов в интерфазе и разделении хромосом при митозе.
Подавились и некоторые гены, нужные при митозе, так как он реже происходит в условиях стресса. Например подавляется ген DNA topoisomerase 2-binding protein 1-A, играющий решающую роль в инициации репликации ДНК и регуляции клеточного цикла. Также снижается экспрессия DNA mismatch repair proteins mutS family domain-containing protein, компонент системы репарации несоответствия ДНК после репликации, так как в условиях стресса реже происходит митоз и репликация. Разные транскрипты одного гена GLYMA_07G261400 (лактоилглутатион-лиаза) по-разному регулируются: в одном случае повысились, в другом понизились. Известно, что белки, содержащие домен лактоилглутатион-лиазы по-разному реагируют на засоление, засуху . Также как и при повышении, среди генов с пониженной регуляцией есть гены, реагирующие на стресс. Например, ген белка 14-3-3. Белки 14-3-3 растений участвуют в сигнальных путях, регулирующих развитие растений и защищающих их от стресса. В растениях комплексы 14–3-3 участвуют в передаче клеточных сигналов, реакции на стресс и регуляции транскрипции , . Понизилась экспрессия гена сульфитоксидазы (SOX); ферменты, связывающие сульфиты, в том числе SOX, по-разному регулируются в условиях абиотического стресса. В частности, они подавляются при солевом стрессе . Понизился ген транскрипционного фактора группы HD-Zip. Считается, что это семейство транскрипционных факторов участвует в регуляции реакции на обезвоживание через разные ветви сигнальной сети, вызывающей обезвоживание, зависимую и независимую от абсцизовой кислоты . Снизилась экспрессия гена белка FAR1-RELATED SEQUENCE. Эти факторы транскрипции участвуют в разных физиологических процессах, таких как передача сигналов фитогормонов, реакция на биотический / абиотический стресс, цветение и старение растений. Воздействие абиотического и биотического стресса может регулировать экспрессию этих генов. Также эти транскрипционные факторы участвуют в передаче световой сигнализации. Участие FRS в обеспечении растений энергией является еще одной возможностью его участия в реакции на стресс, поскольку длительный абиотический стресс подавляет активность фотосистемы II, снижая энергоснабжение растения . Также как и в случае повышения, наиболее сильно понизился ген неохарактеризованного белка, имеющий функцию полиубиквитин-зависимого связывания белков. Это указывает на уменьшение активности системы убиквитин-протеасомной деградации. Снижается экспрессия гена белка, содержащего домен гидроксилазы жирных кислот, участвующего в процессе биосинтеза жирных кислот. Снижается также экспрессия фиолетовой кислотной фосфатазы 36. Физиологические функции фиолетовых кислотных фосфатаз растений разнообразны, включая усвоение фосфора, выработку активных форм кислорода, участие в развитии цветка и биосинтез клеточной стенки .
Снизилась экспрессия субъединицы 1 глутаматдегидрогеназы (GDH). GDH способствует гомеостазу глутамата и играет важную роль на стыке путей ассимиляции углерода и азота. Важно, чтобы растения поддерживали баланс между запасами метаболитов углерода и азота, потому что от этого зависит рост и развитие растений , . Понижение его экспрессии при водном стрессе может указывать на изменение этого баланса, направленное на сохранение энергии и поддержание основных функций клетки в условиях стресса. Подавился ген белка, содержащего домен типа CHY. Эти белки — одни из важных стресс-чувствительных белков растений в реакциях на абиотический стресс посредством деградации белков убиквитинированием . Многие гены этих белков сои могут играть положительную роль регуляторов в реакции на абиотический стресс .
Понизилась экспрессия белка цинкового пальца BRUTUS. Цинковые пальцы выполняют различные биологические функции, такие как связывание с ДНК, РНК, липидами и белками, а также обеспечивают адаптацию к различным биотическим и абиотическим стрессам . Снижается экспрессия RNA polymerase II C-terminal domain phosphatase-like 2 и большой субъединицы фактора инициации транскрипции IIA. Оба участвуют в процессе транскрипции. Снизилась экспрессия белка CHROMATIN REMODELING 4. Он регулирует доступность хроматина, которая необходима как один из механизмов, обеспечивающих активацию генов в нужное время и в нужных типах клеток.
Для валидации списка генов, представленных в таблицах 3 и 4, можно проводить дальнейший анализ методом количественной ПЦР. Лучше всего использовать гены, сильно изменившие экспрессию и в то же время имеющие значения p-value после FDR-коррекции не больше 0,1 (так как они наиболее достоверные и остаются после FDR-коррекции значений p-value). Таким параметрам соответствуют всего два, по одному из каждого списка из 25. Это повысивший экспрессию GLYMA_04G191600 (log2fc≈7,84, p-value 0,081 после коррекции) и понизивший GLYMA_02G208700 (log2fc≈-8,29, p-value 0,006 после коррекции).
Из диаграммы, полученной с помощью ShinyGO (рисунок 2), видно, что наиболее обогащены гены, участвующие в процессе «сборка хиазмы», который обеспечивает правильное разделение гомологичных хромосом во время мейоза , также в репликации ДНК в митохондриях и женском мейотическом клеточном цикле. В целом среди выявленных терминов GO много терминов, относящихся к процессам полового размножения. Это необычный профиль экспрессии для листьев. Это может быть связано с тем, что некоторые гены, участвующие в половом размножении, могут играть двойственную роль в регуляции клеточного цикла, развитии тканей и адаптации к стрессу. Так как мы получили противоречивые данные, то нельзя исключить возможность влияния систематических погрешностей на этапе отбора образцов или подготовки библиотек для секвенирования, что могло повлиять на конечные результаты. Среди других процессов есть такие, которые связаны с убиквитинированием белков, транскрипцией и метаболизмом ксилоглюкана, который присутствует в первичной клеточной стенке высших растений. В других работах среди обогащенных терминов GO также встречаются гены, кодирующие модифицирующие ксилоглюкан ферменты. Это согласуется с наблюдением, что при действии засухи синтезируется больше ксилоглюкана для укрепления клеточной стенки за счет уплотнения и поперечных связей . Наибольшее число генов участвует в процессах клеточного цикла, деления клеточного ядра и органелл, что отражает размер кружка на диаграмме. Вместе с активацией защитных механизмов в ответ на стресс гены, участвующие в росте и развитии, должны быть подавлены .
Среди результатов обогащения GO, полученных с помощью AgriGO, нет терминов, связанных с половым размножением, но много терминов, относящихся к метаболизму углеводов и аминокислот. Это показано в работе , где также в анализе GO большая часть метаболических процессов была связана с метаболизмом углеводов, которые могли обеспечивать энергию, необходимую для процессов при обезвоживании. Процессы метаболизма аминокислот семейства аспартата также способствуют высвобождению энергии, а гены, связанные с метаболизмом аминокислот семейства глутамина и пролина, необходимы для реакции на осмотический стресс .
Несмотря на одинаковый входной набор данных, результаты GO-обогащения, полученных инструментами ShinyGO и AgriGO, сильно различаются. Общий результат — выявление процессов, связанных с метаболизмом ксилоглюкана. Причина расхождений может быть в том, что каждый инструмент использует свои уникальные алгоритмы и базы данных для выявления значимых терминов. Но одним из ключевых факторов, объясняющих несоответствия в результатах GO-обогащения, является различие в уровне покрытия генов, присутствующих в задаваемом списке дифференциально экспрессируемых генов. ShinyGO успешно сопоставил 96% идентификаторов с генами Glycine max, тогда как AgriGO — только около половины. Поэтому результаты GO-обогащения ShinyGO получаются более разнообразными, чем результаты AgriGO, предоставившего более узкий перечень терминов из-за ограниченности сопоставления идентификаторов. Использование комбинированного подхода, включающего оба инструмента, позволяет получать наиболее полное представление о биологических процессах, участвующих в ответе сои на стресс засухи.
В других работах, например, в исследовании реакции сои на засуху сообщалось, что среди генов со сниженной экспрессией наиболее обогащенными терминами GO были активность никотинаминсинтазы, активность биосинтеза никотинаминов и активность фосфатазы. Эти гены вовлечены в биосинтез никотинамина, метаболизм цитокининов и реакцию на окислительный стресс. Однако авторы не выявили никаких существенных терминов обогащения для генов, повышающих экспрессию, так как они в основном были не характеризованы. Другое исследование , направленное на изучение реакции риса на солевой стресс, показало обогащение терминов GO, связанных с общей реакцией на стресс, такие как «ответ на стресс» и «ответ на стимул», а также связанных с транскрипцией генов. В нашем исследовании также обнаружены термины GO, ассоциированные с транскрипцией. Это указывает на важную роль транскрипционных факторов в перестройке профиля транскриптома при реакциях на абиотический стресс . В другом исследовании, посвящённом воздействию соли на растение арабидопсиса , термины GO оказалась ассоциированы с общими реакциями на стрессоры. Наиболее обогащенный термин среди генов с повышенной экспрессией — «реакция на солевой стресс». Вместе с тем, среди генов с пониженной экспрессией выявлены светочувствительные гены, аннотированные как «ответ на интенсивность света» и связанные с восприятием интенсивности освещения. Наконец, в работе по изучению реакции различных видов растений на дефицит воды обнаружены угнетённые процессы, такие как биосинтез воска, формирование клеточной стенки, образование жирных кислот и фосфорилирование белков. Отмечалось усиление таких процессов, как реакция на осмотический стресс, реакция на абсцизовую кислоту, реакция на водную депривацию, активируемый абсцизовой кислотой сигнальный путь, реакция на солевой стресс и реакция на перекись водорода. Эти данные показывают важность общих механизмов приспособления растений к засухе.
5. Заключение
Проведена обработка данных РНК-секвенирования, выявлены гены-маркеры стресса в условиях засухи и засоления: 25 наиболее повысивших экспрессию и 25 наиболее сильно понизивших. Чтобы определить процессы, участвующие в реакции на стресс, сделано GO-обогащение. Анализ GO-обогащения проводился с использованием двух ресурсов: ShinyGO и AgriGO. Полученные результаты показали значительные различия в представленных терминах GO. Из выявленных процессов на ShinyGO многие отвечают за половое размножение, а также за клеточный цикл, среди них также есть процессы, связанные с убиквитинированием белков, транскрипцией и метаболизмом клеточной стенки. Термины GO, связанные с процессом полового размножения необычны для листьев, и причиной такого результата может быть невозможность полностью исключить влияние внешних факторов или потенциальных технических ошибок на этапе подготовки образцов, даже с тщательным соблюдением протоколов сбора материала. Поэтому проводился второй анализ с помощью AgriGO, где получено преобладание процессов, связанных с метаболизмом углеводов и аминокислот. Данные различия результатов GO-обогащения могут быть вызваны различными источниками аннотаций, алгоритмами статистического анализа и способами предварительной фильтрации генов, а также различиями в покрытии идентификаторов. Несмотря на очевидные расхождения, есть частичное совпадение обогащенных терминов, связанных с метаболизмом ксилоглюкана. Для валидации полученного из биоинформатического анализа списка генов, наиболее сильно изменивших регуляцию (маркеров стресса), можно проводить дальнейший анализ методом ПЦР. Кроме того, можно проверить специфичность регуляции этих генов-маркеров на стресс засухи — будут ли они изменять экспрессию при действии других видов стресса.