В данном исследовании был проведен повторный анализ обезличенных, общедоступных данных; набор новых участников, вмешательства или сбор образцов не проводились. В связи с этим для настоящего вторичного анализа не требовалось одобрение институционального комитета по этике и повторного информированного согласия. Сведения об одобрении комитета по этике и информированном согласии для оригинальных исследований были представлены соответствующими поставщиками данных.
Дизайн исследования и логика рабочего процесса
Мы провели ретроспективное биоинформатическое исследование на основе общедоступных данных с использованием наборов данных, размещенных в Gene Expression Omnibus20,21. Новые образцы пациентов, эксперименты на клеточных культурах, модели на животных или лабораторная валидация не проводились. В рабочем процессе использовался дизайн с предварительным определением «якоря» (anchor-first design). Набор GSE17827 применялся для оценки величины эффекта, построения модулей, анализа согласованности между компартментами и анализа градиента симптомов. Независимые наборы данных вводились только после фиксации состава модулей. Рабочий процесс включал анкоринг парных компартментов, картирование по HGNC и фильтрацию белок-кодирующих генов, раздельное построение модулей для слизистой носа и крови, валидацию на соответствующих тканях и в лонгитюдном анализе, клинический бенчмаркинг и анализ робастности. Подтверждающие анализы включали фиксированные тесты на соответствующих тканях и лонгитюдные тесты. Анализы перекрывающихся генов, градиента симптомов, маркерных программ и дисперсии носили исследовательский или описательный характер.
Якорная когорта и основной контраст
GSE17827 содержит профили экспрессии из мазков из носа и цельной крови 26 детей: 9 случаев симптоматического заражения пикорнавирусом, 5 случаев бессимптомного обнаружения пикорнавируса, 6 случаев симптоматического заражения RSV и 6 бессимптомных вирусотрицательных контролей17. По двум случаям RSV отсутствовали образцы крови. Таким образом, полный якорный дизайн включал 50 образцов, а основной контраст «инфицированные против контролей» включал 40 образцов после исключения случаев бессимптомного обнаружения пикорнавируса из построения модулей. Для основного контраста образцы с симптоматическим течением пикорнавирусной и RSV-инфекций были помечены как «инфицированные», а бессимптомные вирусотрицательные образцы — как «контроли». Бессимптомные случаи обнаружения пикорнавируса не рассматривались в качестве контролей, поскольку обнаружение вируса без симптомов биологически отличается от здоровья при отсутствии вируса. Эти образцы были исключены из построения модулей и впоследствии использованы для анализа градиента симптомов.
Картирование и предварительная обработка зондов
Профилирование GSE17827 проводилось на GPL23126. Идентификаторы транскриптных кластеров были сопоставлены с символами генов с использованием файла аннотации Clariom D Human na36, hg38, распространяемого через платформу GEO GPL24539. Были сохранены только символы, одобренные HGNC2. Контрольные зонды, зонды ERCC, некартированные транскриптные кластеры и неодобренные символы были удалены. Несколько транскриптных кластеров, картированных на один и тот же одобренный символ, были объединены по медианному значению экспрессии. После фильтрации и объединения для якорного анализа стало доступно 27 685 генов. Преобразование log₂(x + 1) применялось только в тех случаях, когда 95-й процентиль матрицы экспрессии превышал 50, что указывало на нелогарифмированную шкалу интенсивности. Гены с более чем 20% пропущенных значений были удалены; остальные пропущенные значения были заменены медианой по данному гену. Затем каждый ген был стандартизирован по всем образцам в этом наборе данных как z = (x − mean)/sample standard deviation (ddof = 1). Генам с нулевой дисперсией было присвоено стандартизированное значение 0. Построение модулей было ограничено записями, классифицируемыми как белок-кодирующие гены в полном наборе HGNC.
Анализ величины эффекта и построение модулей
Назальный и гемодинамический компартменты анализировали раздельно. Симптоматические вирусные образцы сравнивали с вирусоотрицательными контролями с использованием g-статистики Хеджеса г стандартизированные средние разности и двусторонние критерии Уэлча23. Хеджирование г был рассчитан как разность средних значений между инфицированной и контрольной группами, деленная на объединенное стандартное отклонение и умноженная на поправку для малых выборок 1 − 3/(4df − 1), где df = нинфицированный + нконтроль − 2. Тесты Уэлча проводились с учетом неравенства дисперсий. Поправка Бенджамини — Хохберга для контроля частоты ложноположительных результатов (FDR) рассчитывалась для всех исследуемых генов в каждом компартменте.24Поскольку ни один белок-кодирующий ген не соответствовал заранее установленному строгому сочетанию показателей FDR < 0,05 и хеджс г > 0,8 в небольшой опорной когорте; гены с положительным эффектом были ранжированы в первую очередь по возрастанию двустороннего t-критерия Уэлча P значение, затем по убыванию g Хеджеса г, при этом символ гена использовался в качестве окончательного определяющего критерия при равенстве показателей. Топ-50 генов образовали каждый первичный модуль. Число пятьдесят было выбрано априори как умеренный размер модуля, обеспечивающий биологическую широту охвата при одновременном ограничении потерь из-за отсутствующих генов на разных платформах. Данный размер модуля не подбирался под внешние значения AUROC. Анализы чувствительности с использованием 10, 25, 50, 10 и 20 генов дали тот же качественный результат разделения якорей. Целью была воспроизводимость на уровне модулей, а не поиск отдельных генов.
Межкомпартментное соответствие
Оценки Hedges g для назальных проб и крови были сопоставлены по генам и сравнены с помощью корреляций Пирсона и Спирмена. Перекрытие между топ-50 модулями в носовых пробах и топ-50 модулями в крови было суммировано по количеству общих генов и индексу Жаккара. Общие гены интерпретировались как предварительные кандидаты на перекрытие рангов между компартментами, а не как валидированные консервативные биомаркеры.
Внешняя валидация на соответствующих тканях
Валидация баллов модулей проводилась с использованием независимых общедоступных наборов данных с интерпретируемыми исходными группами. Валидация по верхним дыхательным путям включала образцы промываний носа при RSV (GSE41374), а также наборы данных по SARS-CoV-2 (GSE152075 и GSE156063). Валидация по крови включала GSE17110 и два отдельно нормализованных блока GSE3890: GPL10558 (36 образцов; 28 RSV и 8 контролей) и GPL6884 (138 образцов; 107 RSV и 31 контроль). Для каждого внешнего набора данных зонды были сопоставлены с символами HGNC, а дублирующиеся символы объединялись по медианному значению экспрессии. В каждом наборе данных применялись те же правила трансформации, фильтрации пропущенных значений, импутации медианы и z-стандартизации по каждому гену, что и для опорного набора данных. Балл модуля рассчитывался как невзвешенное среднее стандартизированных значений экспрессии генов, представленных в данном модуле. AUROC, средняя точность и FDR по тесту Уэлча были представлены в качестве показателей переносимости, а не как оценки клинической диагностической эффективности.
Обширные клинические эталонные наборы данных и наборы данных ex vivo возмущений
Два набора данных цельной крови использовались в качестве клинических эталонов, а не когорт для первичного поиска. В GSE6390 образцы были распределены на основе метаданных о статусе инфекции по группам: вирусная респираторная инфекция (n = 17), бактериальная респираторная инфекция (n = 73) или неинфекционное заболевание (n = 90); образцы без одной из этих однозначных меток были исключены18. В GSE4012 поля диагнозов были сопоставлены с пневмонией, вызванной гриппом А (n = 39), бактериальной пневмонией без гриппа (n = 61), СРОР без пневмонии (n = 40), здоровым контролем (n = 36) или смешанной бактериальной/гриппозной пневмонией (n = 14)19. Образцы со смешанной бактериальной/гриппозной инфекцией были описаны, но исключены из бинарных эталонных сравнений.
Набор GSE53543 содержит 196 профилей PBMC ex vivo от 98 человек. От каждого человека был получен один образец только со средой и один образец, подвергнутый воздействию риновируса 16 в течение 24 h; идентификаторы субъектов подтвердили наличие 98 полных пар. Основной контрольный показатель AUROC был рассчитан для 98 стимулированных и 98 нестимулированных образцов, так как AUROC является метрикой разделения по рангам. Информация о парах субъектов была сохранена в метаданных и использована при проведении парного анализа чувствительности разностей показателей. Данный эксперимент анализировался отдельно от естественной клинической инфекции и не использовался для подтверждения клинических диагностических утверждений.
Были загружены исходные матрицы серий GEO для GSE6390, GSE4012 и GSE53543. Для GSE63990 использовалась GPL571, для GSE4012 — GPL6947, а для GSE53543 — GPL1058. Зонды были сопоставлены с символами, утвержденными HGNC, а дублирующие сопоставления были объединены по медиане экспрессии. Использовались депонированные нормализованные матрицы. Общее правило 95-го процентиля инициировало преобразование log₂(x + 1) только для матриц в нелогарифмическом масштабе. Данные GSE53543 уже были подвергнуты log₂-преобразованию, рангово-инвариантной нормализации и корректировке по дню обработки первоначальными исследователями, поэтому дополнительное логарифмическое преобразование не применялось. После удаления генов с более чем 20% пропущенных значений и заполнения оставшихся пропусков медианными значениями, каждый ген был подвергнут z-стандартизации по всем образцам внутри своего набора данных. Эталонный слой содержал 470 клинических образцов цельной крови и 196 образцов PBMC ex vivo.
Бенчмаркинг эталонных сигнатур
Базовые модули для назальных проб и крови были протестированы с помощью трех невзвешенных эталонных наборов. Официальный набор Pandya из 3 мРНК был перенесен из Дополнительной таблицы 1 оригинального отчета по классификатору5. Компаратор Andres-Terre включал 33 интерферон-ориентированных гена, отобранных из описанной многовирусной сигнатуры3. Компаратор Hallmark включал базовый поднабор из 3 генов интерферона-альфа, связанный с набором MSigDB HALLMARK_INTERFERON_ALPHA_RESPONSE25,26. Полные списки генов для всех трех эталонных наборов приведены в Дополнительных данных 1. Эти эталонные показатели не воспроизводят оригинальные взвешенные классификаторы. Для каждого набора данных показатель рассчитывался как невзвешенное среднее значение доступных z-показателей для каждого гена; требовалось наличие как минимум трех представленных генов, и количество представленных генов было зафиксировано. Сравнения при бенчмаркинге проводились между вирусной и бактериальной инфекциями, вирусной и неинфекционными состояниями, вирусной и бактериальными или неинфекционными состояниями, вирусной и здоровой/контрольной группами, а также бактериальной и здоровой/контрольной группами, если соответствующие группы были доступны. В качестве метрик бенчмаркинга интерпретировались AUROC и средняя точность (average precision). Значения P по тесту Уэлча были скорректированы с использованием процедуры Бенджамини–Хохберга для всех допустимых комбинаций «набор данных — сравнение — модуль» в таблице бенчмаркинга.
Независимая лонгитюдная валидация
Парные наборы данных острого периода и периода выписки GSE9741 и GSE9742 использовались только для лонгитюдной валидации27. Образцы, помеченные в переданных метаданных как одиночная инфекция RSV (RSVsi) или риновирус (hRV), составили основу основного анализа; коинфекции RSV (RSVco) были сохранены только в дополнительных результатах. Метки «острый период» и «выписка» были извлечены из названий образцов. Образцы сопоставлялись по набору данных, компартменту, вирусной группе и идентификатору субъекта, при этом были сохранены только субъекты с данными за обе временные точки. Основная комбинированная группа состояла из 38 пар RSVsi и 30 пар hRV (68 пар на один компартмент).
Для отбора генов, уточнения модулей или настройки пороговых значений наборы данных GSE9741 или GSE9742 не использовались. Эти наборы данных были зарезервированы для внешней валидации. Идентификаторы зондов были сопоставлены с утвержденными символами HGNC, а дублирующиеся символы были объединены по медианному значению экспрессии. Перед расчетом фиксированных показателей модулей GSE17827 были применены те же правила логарифмического преобразования по 95-му процентилю, фильтрация пропущенных значений, импутация медианы и процедуры внутринаборной z-стандартизации для каждого гена.
Для каждой ткани и модуля показатели острого периода и периода выписки были сопоставлены по испытуемому и группе вируса. Для каждой пары рассчитывалась разность между показателем острого периода и показателем периода выписки. Коэффициент Cohen dz рассчитывался как средняя парная разность, деленная на стандартное отклонение этой выборки. Были представлены результаты двухсторонних парных t-критериев и двухсторонних критериев знаковых рангов Уилкоксона, а также AUROC для разделения показателей острого периода и периода выписки. Значения FDR по методу Бенджамини–Хохберга рассчитывались для всех 20 валидных парных t-критериев в полном лонгитюдном выводе (два набора данных, два источника модулей и пять заранее определенных резюме по группам вирусов).
Анализ градиента симптомов и апостериорный анализ
После определения принадлежности к модулям данные о выявлении бессимптомных пикорнавирусов в GSE17827, которые были исключены из основного анализа, использовались только для анализа градиента симптомов. Порядковый клинический балл составлял 0 для вирусоотрицательных контролей, 1 для случаев бессимптомного выявления пикорнавируса и 2 для симптоматической инфекции. Для оценки порядкового тренда использовался коэффициент корреляции Спирмена, а общие различия между тремя группами оценивались с помощью критерия Краскела–Уоллиса. Для сравнения трех пар групп применялись двусторонние апостериорные тесты Манна–Уитни. Поправка Бенджамини–Хохберга применялась отдельно для каждого компартмента к семейству из трех попарных сравнений.
Функциональный анализ обогащения
Гены модулей из носовых выделений и крови были проанализированы с помощью Enrichr через gseapy с использованием библиотек MSigDB Hallmark 2020, Reactome 202 и GO Biological Process 202325,26,28,29,30,31. В Enrichr использовался стандартный метод анализа избыточного представления (over-representation analysis) на основе точного критерия Фишера. Значимыми считались скорректированные по методу Бенджамини — Хохберга P-значения < 0.05, полученные из библиотек. Восемь наиболее значимых терминов для каждого модуля по трем запрашиваемым библиотекам были ранжированы по скорректированному P-значению. Результаты обогащения использовались только для интерпретации.
Анализ робастности, программности маркеров и дисперсии
Анализ робастности проводился для проверки того, зависели ли результаты от размера модуля или случайного выбора. Оценивались модули размером 10, 25, 50, 10 и 20 генов. Для нулевого анализа со случайными генами совокупность подходящих объектов включала белок-кодирующие гены HGNC, представленные в обработанной матрице GSE17827. С помощью генератора случайных чисел NumPy с seed 2026062 было отобрано пятьсот наборов по 50 генов без возвращения. Бутстреп-перевыборка была ограничена 1 0 белок-кодирующими генами с положительным эффектом, имевшими наивысший ранг в исходном якорном анализе для каждого компартмента. В каждой повторной выборке белок-кодирующие гены с положительным эффектом ранжировали по возрастанию двустороннего значения P по критерию Уэлча, а затем по убыванию средней разности. Были отобраны 50 лучших генов. Стабильность генов рассчитывали как количество выборов, деленное на 10. Для каждого компартмента были представлены 20 генов с наибольшей частотой выбора.
Баллы программ-маркеров использовались в качестве описательных средств, а не для оценки клеточных фракций. Были оценены шесть курируемых программ: эпителиальная (EPCAM, KRT8, KRT18, KRT19, MUC1, KRT5, KRT14, SCGB1A1, FOXJ1), моноцитарно-макрофагальная (LYZ, LST1, S10A8, S10A9, FCGR3A, FCGR1A, CD14, MS4A7, CTSS), нейтрофильная (S10A8, S10A9, MPO, ELANE, CEACAM8, FCGR3B, OLFM4, MMP8), T/NK (CD3D, CD3E, TRAC, NKG7, GNLY, PRF1, GZMB, KLRD1, IL7R), B/плазматическая (MS4A1, CD79A, CD79B, MZB1, JCHAIN, IGHG1, SDC1) и миелоидный интерфероновый ответ (SIGLEC1, IFI27, IFI4L, ISG15, MX1, OAS1, RSAD2, IFIT3). Балл каждой программы рассчитывался как среднее значение доступных z-показателей для отдельных генов, при этом требовалось наличие как минимум трех представленных генов. Значения P по критерию Спирмена были скорректированы с использованием процедуры Бенджамини–Хохберга для всего семейства корреляций «набор данных–модуль–программа». Однофакторный показатель eta-квадрат рассчитывался как сумма квадратов между группами, деленная на общую сумму квадратов для каждой пары «модуль–фактор». Строки, в которых отсутствовал соответствующий фактор, исключались из этого расчета. Данный анализ носил описательный характер и не учитывал взаимно коррелирующие факторы.
Воспроизводимость
Все виды анализа выполнялись с использованием скриптов на языке Python 3.12.13 с применением программных пакетов и версий, указанных в Таблице материалов. Для рандомизированных процедур использовалось фиксированное начальное число (seed) 2026062. Публичные матрицы экспрессии из GEO обрабатывались с использованием единообразной структуры проекта, в которой были разделены необработанные входные данные, обработанные данные, результаты статистического анализа, таблицы и рисунки. Для обеспечения возможности независимой проверки были сохранены определения модулей, записи о курации образцов, оценки величины эффекта, статистика валидации, результаты обогащения, анализ робастности и исходные данные для рисунков. Код анализа, спецификации зависимостей и вспомогательные производные данные доступны в соответствии с разделом «Доступность данных».