La réanalyse des transcriptomes publics révèle des signatures immunitaires partagées entre le trouble dépressif majeur et la dermatomyosite avec un contexte unicellulaire

Résumé

Cette étude visait à utiliser une réanalyse bioinformatique intégrative de jeux de données publics GEO, combinée à une contextualisation unicellulaire, pour identifier des gènes candidats partagés entre trouble dépressif majeur et dermatomyosite, et caractériser leur répartition à travers des sous-ensembles cellulaires immunitaires dans un ensemble de données unicellulaire lié à la dermatomyosite.

Résumé

Cette étude visait à identifier des signaux transcriptomiques partagés candidats entre trouble dépressif majeur et dermatomyosite par une réanalyse bioinformatique intégrative de jeux de données publiques GEO avec contextualisation unicellulaire. Le flux de travail analytique comprenait l’analyse pondérée du réseau de co-expression génétique (WGCNA) pour l’identification des modules clés, l’ontologie génétique (GO) et les analyses d’enrichissement de la Kyoto Encyclopedia of Genes and Genomes (KEGG) pour la caractérisation fonctionnelle, l’analyse de réseau basée sur GeneMANIA et une plateforme de visualisation de réseau pour la priorisation des gènes candidats, ainsi que l’évaluation de 113 modèles d’apprentissage automatique combinés avec SHapley Additive ExPlanations (SHAP) pour la sélection des caractéristiques diagnostiques. L’analyse d’enrichissement des ensembles de gènes (GSEA), l’analyse d’infiltration immunitaire et la contextualisation basée sur l’ARN-seq de cellules uniques ont ensuite été réalisées afin de mieux caractériser le contexte cellulaire lié au système immunitaire des signaux identifiés. L’intégration de jeux de données GEO liés à la dermatomyosite a permis d’identifier 570 gènes exprimés différemment, dont 33 gènes candidats partagés ont été obtenus via WGCNA. L’enrichissement fonctionnel et les analyses en réseau ont mis en lumière la défense immunitaire, la cytotoxicité et les voies comprenant PPAR, IL-17 et traitement antigénique, avec ELANE, PPBP et CTSG émergeant comme des nœuds fortement connectés. La priorisation des fonctionnalités basée sur l’apprentissage automatique a conservé 8 gènes sélectionnés par le modèle candidat, à savoir KIF4A, OLR1, KIR2DL4, KRT23, KIR3DS1, AZU1, SCG5 et LRRC37E. L’analyse d’infiltration immunitaire a associé ces gènes communs avec les lymphocytes T régulateurs (Tregs), les mastocytes au repos, les cellules dendritiques au repos, ainsi que les macrophages à activation classique (M1) et alternativement activée (M2). La contextualisation de l’ARN-seq unicellulaire a également suggéré que les sous-ensembles de cellules T CD8⁺ présentant différents états de score du gène candidat présentaient des schémas de communication intercellulaires distincts. Parmi celles-ci, l’axe MIF–(CD74+CD44) et les signaux provenant des cellules T à mémoire naïve/centrale étaient des caractéristiques notables nécessitant une validation supplémentaire. Dans l’ensemble, cette étude a identifié des signaux transcriptomiques partagés candidats entre trouble dépressif majeur et dermatomyosite et a mis en lumière des contextes cellulaires liés au système immunitaire qui justifient une validation supplémentaire dans les cohortes véritablement comorbides.

Introduction

La dermatomyosite est une maladie auto-immune systémique chronique caractérisée par une atteinte inflammatoire de la peau et du muscle squelettique, manifestée cliniquement par une faiblesse musculaire proximale symétrique et des lésions cutanées distinctives, et, dans les cas graves, un dysfonctionnement^{multiorgane 1}. L’accumulation de preuves cliniques met en lumière que les patients atteints de dermatomyosite présentent fréquemment des comorbidités psychiatriques, notamment le trouble dépressif^{majeur 2,3,4}. La pathogenèse du trouble dépressif majeur associé à la dermatomyosite est multifactorielle, résultant d’une interaction complexe entre détresse psychosociale, dysrégulation neuroendocrinienne et immunoinflammation systémique. La douleur persistante, la fatigue et une faiblesse musculaire progressive peuvent gravement nuire à la fonction physique et à la qualité de vie. Ces fardeaux peuvent entraîner des perturbations des rôles sociaux, un stress psychologique chronique et une autonomie^{réduite 5}. De plus, une exposition prolongée aux glucocorticoïdes peut perturber l’axe hypothalamo-hypophyse-surrénale (HPA) et altuer la plasticité de l’hippocampe, augmentant ainsi la vulnérabilité aux manifestations liées à^{la dépression 6}. De plus, l’activation immunitaire soutenue et l’inflammation systémique dans la dermatomyosite sont de plus en plus reconnues comme des facteurs potentiels contributeurs aux symptômes liés à^{la dépression 7}. Ces médiateurs périphériques peuvent influencer le système nerveux central en modulant le métabolisme des neurotransmetteurs et la neuroplasticité, liant ainsi l’auto-immunité systémique aux manifestations neuropsychiatriques 7,8,9,10.

Il est important de noter que cette justification n’implique pas que l’une ou l’autre maladie soit biologiquement homogène. La dermatomyosite comprend des sous-ensembles cliniquement et sérologiquement distincts, incluant des phénotypes associés à l’anti-MDA5 et à l’anti-TIF1-γ avec des profils inflammatoires et cliniques^différents 11,12,13. Le trouble dépressif majeur est également de plus en plus reconnu comme une affection hétérogène, et les preuves actuelles soutiennent l’existence d’un sous-type immunitaire-inflammatoire plutôt qu’une seule signature inflammatoire^{universelle 14,15}. En conséquence, la présente étude n’a pas été conçue pour supposer un programme moléculaire partagé universel, mais plutôt pour dépister les signaux transcriptomiques immunoassociés chevauchants des candidats détectables au niveau de la cohorte à travers des ensembles de données publics indépendants.

Au-delà du stress psychosocial et de l’exposition au traitement, un lien biologiquement plus testable entre trouble dépressif majeur et dermatomyosite est la dysrégulation immunitaire-inflammatoire partagée. Le trouble dépressif majeur est une condition hétérogène et ne doit pas être supposée avoir un profil transcriptomique universel unique. Cependant, des preuves convergentes soutiennent un sous-type lié à l’inflammation du trouble dépressif majeur, et des études transcriptomiques périphériques ont identifié une dysrégulation des voies immunitaires innées, liées aux neutrophiles, de l’interféron et du complément chez des sous-groupes d’individus^atteints 14,16,17. Parallèlement, la signalisation de stress liée à la MACK a également été impliquée dans les phénotypes^{dépressifs 18}. La dermatomyosite, en revanche, est une maladie auto-immune bien reconnue à l’interféron, et des études transcriptomiques dans le sang et les tissus affectés ont constamment démontré l’activation de l’interféron de type I et de programmes immunitaires inflammatoires plus larges ; des analyses multi-omiques récentes ont également mis en lumière l’activité des voies ERK et p38 liées à MAPK dans la dermatomyosite 19,20,21. Collectivement, ces résultats fournissent une justification biologiquement plausible pour examiner si un sous-ensemble de signaux transcriptomiques associés au système immunitaire peut recouper le trouble dépressif majeur et la dermatomyosite dans des ensembles de données publics indépendants.

Malgré ces observations, la base moléculaire sous-jacente au chevauchement entre le trouble dépressif majeur et la dermatomyosite reste insuffisamment comprise. Il est important de noter que les ensembles de données publics actuellement disponibles ne fournissent pas une véritable cohorte de patients diagnostiqués simultanément avec un trouble dépressif majeur et une dermatomyosite. Ainsi, plutôt que d’analyser directement la dépression chez les patients atteints de dermatomyosite, la présente étude a été conçue pour identifier les signaux transcriptomiques partagés candidats à travers des ensembles de données publics distincts sur le trouble dépressif majeur et la dermatomyosite via un cadre de réanalyse bioinformatique^{intégrative 22}. Plus précisément, des ensembles de données transcriptomiques accessibles au public ont été analysés à l’aide d’une analyse d’expression différentielle, d’une analyse pondérée des réseaux de co-expression génique (WGCNA), d’une analyse d’enrichissement fonctionnel, d’une analyse basée sur des réseaux et d’une priorisation de caractéristiques basée sur l’apprentissage automatique afin d’identifier les gènes et voies candidats^{inter-maladies 23}. De plus, un ensemble de données unicellulaire lié à la dermatomyosite a été analysé pour contextualiser ces gènes candidats au niveau des cellules immunitaires. Comme illustré à la Figure 1, le flux de travail analytique global est résumé dans un organigramme par étapes. Plutôt que d’établir un mécanisme définitif de comorbidité, cette étude visait à générer un cadre guidé par des hypothèses pour identifier les signaux moléculaires partagés par les candidats entre trouble dépressif majeur et dermatomyosite.

En conséquence, cette étude a adopté un cadre de priorisation par étapes. Les modules de co-expression associés à la maladie ont d’abord été identifiés séparément dans les ensembles de données sur le trouble dépressif majeur et la dermatomyosite, et leur chevauchement a été utilisé pour définir les signaux partagés entre maladies candidates. Ces candidats ont ensuite été contextualisés fonctionnellement par enrichissement et analyse de réseau basée sur GeneMANIA, priorisés dans un cadre de classification centré sur la dermatomyosite utilisant des méthodes d’apprentissage automatique, puis enfin examinés dans un ensemble de données unicellulaire lié à la dermatomyosite pour fournir une contextualisation cellulaire.

Figure 1 : Organigramme du processus de collecte et d’analyse des données. Veuillez cliquer ici pour voir une version agrandie de cette figure.

Plusieurs approches alternatives ont été utilisées pour étudier le chevauchement moléculaire entre les maladies. L’intersection simple des listes de gènes (DEG) exprimées différenciellement est simple sur le plan informatique mais manque des informations contextuelles au niveau des modules fournies par l’analyse de co-expression et est sensible à des seuils arbitraires de changement de pli et de valeurs P. La méta-analyse traditionnelle regroupe les tailles d’effet entre les études d’une même maladie mais n’est pas conçue pour identifier des signaux partagés entre deux affections distinctes. Le flux de travail actuel intègre plusieurs couches analytiques complémentaires — chevauchement de modules de co-expression, enrichissement fonctionnel, analyse de réseau, priorisation des fonctionnalités basée sur l’apprentissage automatique, déconvolution immunitaire et contextualisation à cellule unique — chacune remplissant un rôle distinct dans un cadre de priorisation séquentielle. Cette conception à plusieurs niveaux aide à réduire le nombre de gènes candidats par étapes et offre une contextualisation biologique validée à plusieurs niveaux. Le protocole s’applique à toute paire de maladies pour lesquelles il existe des ensembles de données transcriptomiques en vrac et, éventuellement, monocellulaires, disponibles publiquement, en particulier lorsque les cohortes véritables comorbides ne sont pas disponibles. Cependant, le flux de travail est observatoire et n’intègre pas de cadres formels d’inférence causale ; Tous les résultats doivent être interprétés comme générateurs d’hypothèses et nécessitent une validation expérimentale indépendante.

Dans l’ensemble, le flux de travail analytique a été conçu comme une stratégie de priorisation séquentielle plutôt que comme un cadre d’inférence causale directe. Chaque étape avait un objectif distinct : chevauchement de modules basés sur WGCNA pour l’identification de signaux partagés candidats, analyse d’enrichissement/réseau pour la contextualisation biologique, apprentissage automatique pour la priorisation des caractéristiques dans la classification liée à la dermatomyosite, et analyse unicellulaire pour la contextualisation au niveau des types cellulaires.

Protocole

Cette étude n’a utilisé que des ensembles de données publiques et désidentifiées provenant de la base de données Gene Expression Omnibus (GEO). Comme le travail impliquait une analyse secondaire de données publiques existantes et n’incluait pas de contact direct des participants, d’intervention ou d’accès à des informations personnelles identifiables, une approbation supplémentaire du comité d’éthique ni un consentement éclairé n’étaient pas requis.

Sources de données et prétraitement

Tous les ensembles de données sur l’expression génique et les jeux de données unicellulaires ont été obtenus à partir de la base de données^{GEO 24}. Pour le trouble dépressif majeur, GSE98793 de données a été utilisé, comprenant des échantillons de sang périphérique de 128 patients et 64 témoins sains. Pour la dermatomyosite, les ensembles de données ont été sélectionnés selon des critères prédéfinis, incluant le profilage de l’expression de Homo sapiens, les groupes de maladies et de témoins clairement identifiables, l’annotation disponible de plateforme pour la cartographie sonde-gène, et l’adéquation à l’analyse de découverte ou de validation. Lorsqu’une série GEO contenait plusieurs sous-types de myopathie inflammatoire, seuls des échantillons de dermatomyosite et de témoins normaux ont été prélevés pour la présente étude. GSE1551, GSE46239 et GSE128470 ont été utilisés comme ensembles de données de découverte/entraînement, tandis que GSE5370, GSE39454 et GSE11971 ont été utilisés comme ensembles de données de validation indépendants. Les ensembles de données sur la dermatomyosite analysés dans cette étude provenaient principalement de tissus musculaires ou cutanés affectés plutôt que du sang périphérique. Les données unicellulaires pour la dermatomyosite proviennent de l’ensemble de données GSE190510.

Les matrices d’expression brutes étaient téléchargées depuis la base de données GEO avec les fichiers d’annotation de la plateforme correspondants. Les identifiants des sondes étaient mappés aux symboles officiels des gènes selon l’annotation GPL fournie par le fabricant. Les sondes qui ne pouvaient être assignées sans ambiguïté à un seul symbole génique officiel ont été retirées. Lorsque plusieurs sondes étaient mappées sur le même gène, elles étaient effondrées au niveau du gène en utilisant la valeur d’expression moyenne implémentée par la fonction « avereps » dans le package limma, générant ainsi une matrice d’expression gène par échantillon.

Pour réduire le biais dépendant de l’intensité et stabiliser la variance, la transformation log2 a été appliquée lorsque approprié selon la distribution des valeurs d’expression. La normalisation entre tableaux était alors effectuée à l’aide de la fonction « normalizeBetweenArrays » dans le package limma. Les valeurs manquantes, lorsqu’elles étaient présentes, étaient imputées en utilisant l’imputation de K-plus proche voisin. Pour les ensembles de données intégrés d’entraînement à la dermatomyosite, la correction par lots a été réalisée à l’aide de la fonction « ComBat » dans le package SVA, l’origine du jeu de données/plateforme étant traitée comme variable batch et groupe échantillon (dermatomyosite versus contrôle sain) inclus dans la matrice de conception afin de préserver la variation biologique d’intérêt lors de l’ajustement des lots.

Toutes les analyses étaient réalisées en R à l’aide d’un environnement de développement intégré pour R sur un système d’exploitation de bureau. Le paquet limma était utilisé pour le résumé et la normalisation des sondes. Le package sva était utilisé pour la correction batch de ComBat. Les valeurs manquantes ont été imputées en utilisant l’imputation de K-plus proche voisin avec k = 10.

Analyse pondérée du réseau de co-expression génique

L’analyse pondérée du réseau de co-expression génique (WGCNA) a été réalisée séparément pour les ensembles de données sur le trouble dépressif majeur et la dermatomyosite à l’aide du package WGCNA^{R 25,26}. Les échantillons étaient regroupés hiérarchiquement à l’aide de flashClust pour identifier les valeurs aberrantes ; Les échantillons dépassant une hauteur de dendrogramme de 100 et les gènes situés dans les 25 % les plus bas de la variance ont été exclus. Pour chaque réseau, une puissance de seuil doux (β) était sélectionnée à l’aide de pickSoftThreshold afin d’obtenir une topologie approximativement sans échelle (R² > 0,8). La matrice d’adjacence a été transformée en matrice de chevauchement topologique (TOM), et les modules ont été identifiés via la coupe dynamique de l’arbre avec une taille minimale de module de 60 et une hauteur de coupe de fusion de 0,25²⁷. Le paquet WGCNA R était utilisé avec flashClust pour le clustering hiérarchique. La graine aléatoire a été réglée à 12345 pour la reproductibilité. Les eigengènes du module ont été corrélés au statut pathologique grâce à la corrélation de Pearson, avec des valeurs P ajustées par la méthode de Benjamini–Hochberg. Pour chaque maladie, le module montrant l’association la plus forte et la plus significative avec le statut de maladie a été conservé comme module clé associé à la maladie. Le chevauchement entre les gènes clés du module de données sur le trouble dépressif majeur et ceux du jeu de données sur la dermatomyosite a été défini comme le jeu de gènes candidats partagés pour les analyses en aval. Une analyse d’expression différentielle de la cohorte intégrée de dermatomyosite a été réalisée séparément afin de caractériser les modifications transcriptionnelles liées à la dermatomyosite.

Analyse de l’enrichissement fonctionnel

L’analyse d’enrichissement par ontologie génique (GO) a été réalisée à l’aide de R. Les symboles géniques ont été convertis en identifiants Entrez via org. Hs.eg.db, et des termes GO significativement enrichis (p < 0,05) ont été identifiés en utilisant enrichGO dans clusterProfiler. Pour la visualisation multidimensionnelle des résultats, des diagrammes en barres et des diagrammes à bulles ont été générés à l’aide du package enrichplot, tandis qu’un diagramme circulaire a été construit avec le package circlize pour afficher les catégories GO, le nombre de gènes et les facteurs d’enrichissement. Les légendes ont été ajoutées avec le package ComplexHeatmap. L’analyse d’enrichissement des voies de gènes différenciellement exprimés par l’Encyclopédie des gènes de Kyoto Encyclopedia of Genes and Genomes (KEGG) a également été réalisée dans R. Les symboles génétiques ont été convertis en identifiants Entrez basés sur l’organisation. Hs.eg.db base de données et des voies significativement enrichies (FDR < 0,05) ont été identifiées à l’aide de la fonction enrichKEGG du package^{clusterProfiler} 28,29,30,31. Les résultats d’enrichissement ont été visualisés à l’aide de diagrammes en barres et de bulles.

Analyse de réseau d’association fonctionnelle basée sur GeneMANIA

Sur la base des gènes partagés précédemment identifiés, un réseau d’association fonctionnelle basé sur GeneMANIA a été construit pour explorer le contexte d’interaction entre ces gènes et leurs partenaires apparentés. La liste des gènes a été soumise à GeneMANIA en utilisant Homo sapiens comme espèce de référence. GeneMANIA intègre plusieurs types de preuves, notamment la co-expression, les interactions physiques, les voies, la co-localisation, les interactions génétiques et les domaines protéiques partagés. Le réseau résultant a été exporté et importé dans une plateforme de visualisation réseau pour la visualisation et l’analyse. L’analyse topologique du réseau a ensuite été réalisée sur une plateforme de visualisation réseau pour la visualisation et l’analyse afin d’identifier les nœuds candidats hautement^{connectés 32, 33, 34}.

Construction de modèles de diagnostic basée sur l’apprentissage automatique

Plusieurs algorithmes d’apprentissage automatique ont été utilisés pour la classification diagnostique, notamment Random Forest (RF), Support Vector Machine (SVM), Linear Discriminatant Analysis (LDA), Naive Bayes, Gradient Boosting Machine (GBM), XGBoost, glmBoost, Elastic Net (Enet), Ridge, Least Absolute Shrinkage and Selection Operator (LASSO), Stepglm généralisé Linear Model (Stepglm) et Partial Least Least Squares Generalized Linear Model (plsRglm⁾³⁵. Un cadre de modélisation à deux étapes a été appliqué pour générer 113 combinaisons de modèles candidats. Lors de la première étape, l’algorithme initial a été utilisé pour le criblage variable dans la cohorte de formation ; À la deuxième étape, les variables conservées ont été utilisées pour ajuster un modèle de classification diagnostique. Les modèles avec ≤5 variables sélectionnées ont été exclus de la comparaison ultérieure. Les ensembles de données combinés sur la dermatomyosite ont servi de cohorte d’entraînement, avec des étiquettes définies comme dermatomyosite versus témoins sains, tandis que les cohortes de validation indépendantes étaient utilisées pour l’évaluation externe de la performance. Le rééchantillonnage interne et l’ajustement étaient spécifiques à l’algorithme : les modèles basés sur glmnet (LASSO, Ridge et Elastic Net) utilisaient une validation croisée 10 fois pour sélectionner lambda.min ; GBM utilisait une validation croisée interne à 10 fois pour déterminer le nombre optimal d’arbres ; XGBoost a utilisé un rééchantillonnage en 5 fois pour sélectionner le dernier tour de boost selon la perte logarithmique minimale du test ; glmBoost utilisait une validation interne interne basée sur cvrisk pour déterminer l’itération d’arrêt ; et LDA était intégrée au cadre de validation croisée CARET. Pour les algorithmes sans étapes d’ajustement explicites dans la mise en œuvre actuelle, des paramètres fixes ou par défaut étaient utilisés. Pour réduire la fuite d’informations, la sélection de caractéristiques, l’ajustement des modèles et l’ajustement interne ont été réalisés uniquement en utilisant la cohorte de formation, tandis que la ou les cohortes de validation ont été utilisées uniquement pour la prédiction indépendante et l’évaluation de performance basée sur l’AUC. Le package caret était utilisé pour la gestion des flux de travail en apprentissage automatique, avec glmnet, randomForest, e1071, gbm, xgboost, mboost, plsRglm et MASS pour des algorithmes individuels. L’analyse SHAP a été réalisée à l’aide du package shapviz. La graine aléatoire était réglée à 12345 avant chaque montage du modèle. Les modèles avec moins de 5 caractéristiques sélectionnées étaient exclus. L’interprétabilité du modèle et la contribution au niveau des gènes ont été ensuite évaluées à l’aide des explications additives SHapley (SHAP), et les gènes les plus informatifs ont été priorisés comme caractéristiques sélectionnées par le modèle candidat pour l’interprétation biologique en aval.

Évaluation de la performance diagnostique

Les courbes de caractéristiques de fonctionnement du récepteur (ROC) ont été générées à l’aide du package R « pROC » pour évaluer la performance diagnostique des biomarqueurs candidats. Les niveaux d’expression et la précision prédictive des marqueurs candidats ont été validés dans des ensembles de données indépendants (GSE5370, GSE11971 et GSE39454). La performance du modèle a été ensuite évaluée à l’aide de matrices de confusion. L’expression différentielle des gènes des modules clés a été visualisée à l’aide de volcans et de boîtes de tableaux, et des courbes ROC ont été construites pour évaluer la valeur diagnostique de chaque gène.

Analyse d’enrichissement d’ensembles géniques

Pour explorer les changements fonctionnels coordonnés associés aux signaux transcriptomiques partagés candidats, l’analyse d’enrichissement des ensembles de gènes (GSEA) a été réalisée à l’aide de clusterProfiler^36,37. Les données d’expression génique issues de la dermatomyosite et des échantillons témoins ont été classées selon l’expression différentielle. Des ensembles de gènes prédéfinis correspondant aux voies KEGG (c2.cp.kegg.Hs.symbols.gmt) ont été utilisés pour évaluer si les gènes de chaque voie présentaient une tendance coordonnée à la hausse ou à la baisse. La signification statistique a été définie comme P < 0,05.

Analyse de l’infiltration des cellules immunitaires

La matrice de dermatomyosite normalisée, transformée en log2 et corrigée en lot a été utilisée pour la déconvolution immunitaire. L’algorithme CIBERSORT a été appliqué pour estimer l’abondance relative des sous-types de cellules immunitaires à l’aide de la matrice de référence^{LM22 38}. Des échantillons avec la déconvolution P < 0,05 ont été conservés pour une analyse en aval. Les différences dans les proportions inférées de cellules immunitaires entre les groupes ont été visualisées à l’aide de boîtes graphiques, et une analyse de corrélation de Spearman a été réalisée pour évaluer les associations entre les sous-ensembles de cellules immunitaires et les gènes candidats partagés.

Analyse du séquençage de l’ARN unicellulaire pour la contextualisation cellulaire

Des analyses RNA-seq en cellules uniques ont été réalisées en R à l’aide de Seurat. Harmony a été utilisé pour la correction par lots, DoubletFinder pour la détection de doublets, celda/decontX pour l’estimation de l’ARN ambiant, Monocle pour l’analyse de trajectoire pseudo-temps, CellChat pour l’analyse de la communication cellulaire-cellule, AUCell pour l’évaluation de l’activité des ensembles génétiques, et GSVA pour le score ssGSEA. Les matrices de comptage brut ont été importées dans des objets Seurat avec les paramètres min.cellules = 5 et min.caractéristiques = 300. Des indicateurs de contrôle qualité, incluant les proportions des gènes mitochondriaux, ribosomiaux et d’hémoglobine, ont été calculés pour chaque cellule. Les cellules n’étaient conservées que si elles remplissaient tous les critères suivants : nFeature_RNA > 500, nCount_RNA < 5 000, percent_mito < 25, percent_ribo > 3 et percent_hb < 1. Les gènes détectés dans moins de 3 cellules ont été exclus. De plus, MALAT1 et les gènes mitochondriaux ont été retirés avant l’analyse en aval. Après le filtrage initial, des doublets ont été identifiés dans chaque échantillon à l’aide de DoubletFinder, avec des PC = 1:30 et pN = 0,25 ; Les taux de doublet attendus ont été fixés selon le nombre de cellules spécifiques à chaque échantillon (<4 000 cellules : 2,5 % ; 4 000–8 000 cellules : 5 % ; >8 000 cellules : 6,5 %). Seuls les débardeurs furent conservés. La contamination à l’ARN ambiant a été estimée en outre à l’aide de decontX, et des cellules avec des scores de contamination < 0,2 ont été conservées.

Les données filtrées ont été normalisées à l’aide de la méthode LogNormalize avec un facteur d’échelle de 10 000, suivies d’identification des gènes variables, d’une mise à l’échelle des données et d’une analyse des composantes principales. Les effets batch entre échantillons ont été corrigés en utilisant Harmony avec orig.ident comme variable batch. Les 15 premières dimensions Harmony ont été utilisées pour la visualisation UMAP et la construction de graphes voisins. Le clustering a été réalisé à l’aide de FindNeighbors et FindClusters, et le résultat final du clustering a été défini à une résolution de 0,05. Les types cellulaires ont été annotés manuellement selon les gènes marqueurs canoniques ainsi que les résultats FindAllMarkers³⁹.

Pour la contextualisation fonctionnelle en aval, l’activité des gènes candidats a été évaluée au niveau de la cellule unique, et le sous-ensemble pertinent des cellules immunitaires a été soumis à des analyses de trajectoire et de communication intercellulaire. L’analyse du pseudotemps a été réalisée à l’aide de Monocle avec réduction de dimensionnalité basée sur DDRTree, suivie d’un ordre cellulaire. L’analyse de la communication cellulaire-cellule a été réalisée à l’aide de CellChat avec la base de données ligand-récepteur humain, limitée à la catégorie de signalisation sécrétée, et les communications impliquant moins de 10 cellules ont été filtrées.

Pour chaque cellule, l’activité du gène candidat a été quantifiée à l’aide de trois approches complémentaires : AUCell, ssGSEA et AddModuleScore. Les scores AUCell ont été calculés à partir de matrices de classement génétique, et les scores ssGSEA ont été générés à l’aide du cadre GSVA. AddModuleScore a été calculé à l’aide de la fonction intégrée Seurat. Les valeurs résultantes d’AUCell, ssGSEA et AddModuleScore ont ensuite été combinées en une seule matrice de scores. Chaque type de score a d’abord été standardisé par transformation Z-score, puis redimensionné à une plage de 0 à 1 en utilisant la normalisation min–max. Le score composite final (« Score ») pour chaque cellule a été défini comme la somme des trois scores normalisés :

Score = AUCell normalisé + ssGSEA normalisé + AddModuleScore normalisé.

Pour les analyses de sous-groupes en aval, le sous-ensemble des cellules T CD8⁺ a été extrait, et les cellules ont été dichotomisées selon la valeur médiane de Scoring au sein de ce sous-ensemble. Les cellules avec des valeurs de score supérieures à la médiane ont été attribuées au groupe High_Hub_genes, tandis que les cellules restantes ont été assignées au groupe Low_Hub_genes.

Résultats

Identification des gènes candidats partagés entre le trouble dépressif majeur et la dermatomyosite

Après la fusion des données, la normalisation et la correction par lots de jeux de données GEO liés à la dermatomyosite (Figure 2A,B), un total de 570 gènes exprimés différemment ont été identifiés (Figure 2C,D), comprenant 517 gènes à la hausse et 53 à la baisse.

Discussion

La dermatomyosite est une maladie auto-immune systémique chronique avec une atteinte importante de la peau et des muscles, et les observations cliniques accumulées suggèrent que les patients atteints de dermatomyosite peuvent également ressentir une charge psychiatrique importante, y compris des symptômes compatibles avec un trouble dépressif majeur. Dans ce contexte, la présente étude a appliqué un cadre de réanalyse bioinformatique intégrative pour identifier les signaux transcriptomiq...

Matériaux

Liste des matériaux utilisés dans cet article

Nom	Entreprise	Numéro de catalogue	Commentaires
AddModuleScore	Fonction Seurat	version 4.4.0	Calcul du score du module au sein de Seurat RRID: NA
AUCell	Bioconductor	version 1.32.0	Scorage de l'activité des ensembles de gènes à cellule unique RRID: SCR_021327
caret	CRAN	version 7.0.1	Support du workflow d'apprentissage automatique RRID: SCR_022524
celda / decontX	Bioconductor	version 1.24.0	Estimation de la contamination par l'ARN ambiant RRID: NA
CellChat	GitHub / CellChat	version 2.2.0	Analyse de la communication cellule-cellule RRID: SCR_021946
CIBERSORT / LM22 signature matrix	CIBERSORT	LM22	Estimation de l'infiltration des cellules immunitaires RRID: NA
clusterProfiler	Bioconductor	version 4.12.6	Analyse d'enrichissement fonctionnel RRID: SCR_016884
Cytoscape	Cytoscape Consortium	version 3.10	une plateforme de visualisation de réseau pour la visualisation et l'analyse RRID: SCR_003032
DoubletFinder	GitHub / McGinnis Lab	version 2.0.4	Détection des doublets dans les ensembles de données à cellule unique RRID: NA
e1071	CRAN	version 1.7.16	Modélisation de la machine à vecteurs de support et de Bayes naïf RRID: NA
gbm	CRAN	version 2.2.2	Modélisation de Gradient Boosting Machine RRID: NA
Gene Expression Omnibus (GEO) database	National Center for Biotechnology Information (NCBI)	GSE98793	Ensemble de données du transcriptome de trouble dépressif majeur RRID: NA
Gene Expression Omnibus (GEO) database	NCBI	GSE1551	Ensemble de données d'apprentissage de dermatomyosite; échantillons de biopsie musculaire squelettique RRID: NA
Gene Expression Omnibus (GEO) database	NCBI	GSE46239	Ensemble de données d'apprentissage de dermatomyosite; échantillons de biopsie cutanée RRID: NA
Gene Expression Omnibus (GEO) database	NCBI	GSE128470	Ensemble de données d'apprentissage de dermatomyosite; échantillons de dermatomyosite extraits du groupe de myopathies inflammatoires RRID: NA
Gene Expression Omnibus (GEO) database	NCBI	GSE5370	Ensemble de données de validation indépendant de dermatomyosite; échantillons de muscle adulte non traités RRID: NA
Gene Expression Omnibus (GEO) database	NCBI	GSE11971	Ensemble de données de validation indépendant de dermatomyosite RRID: NA
Gene Expression Omnibus (GEO) database	NCBI	GSE39454	Ensemble de données de validation indépendant de dermatomyosite; échantillons de dermatomyosite extraits du groupe de myopathies inflammatoires RRID: NA
Gene Expression Omnibus (GEO) database	NCBI	GSE190510	Ensemble de données de RNA-seq à cellule unique liées à la dermatomyosite RRID: NA
GeneMANIA	University of Toronto / GeneMANIA	version du serveur web accédée dans cette étude	Construction de réseau d'association fonctionnelle RRID: RRID:SCR_005709
glmnet	CRAN	version 4.1.8	Modélisation LASSO, Ridge et Elastic Net RRID: NA
GSVA	Bioconductor	version 2.0.7	Scorage ssGSEA RRID: NA
Harmony	CRAN	version 1.2.4	Correction des lots pour l'intégration des données à cellule unique RRID: NA
limma	Bioconductor	version 3.60.6	Analyse de l'expression différentielle, résumé des sondes et utilitaires de normalisation RRID: SCR_010943
MASS	CRAN	version 7.3.61	Analyse discriminatoire linéaire RRID: NA
mboost	CRAN	version 2.9.11	Modélisation glmBoost RRID: NA
Monocle	Bioconductor	version 2.38.0	Analyse de trajectoire pseudo-temps RRID: SCR_016339
org.Hs.eg.db	Bioconductor	version 3.19.1	Base de données d'annotation des gènes humains RRID: NA
plsRglm	CRAN	version 1.5.1	Modélisation par régression PLS généralisée RRID: NA
pROC	CRAN	version 1.18.5	Analyse de la courbe ROC RRID: SCR_024286
R statistical software	R Foundation for Statistical Computing	version 4.4.2	Environnement principal de calcul statistique RRID: SCR_001905
randomForest	CRAN	version 4.7.1.2	Modélisation Random Forest RRID: SCR_015718
RStudio	Posit Software, PBC	version 2024.4.1.748	Environnement de développement intégré pour R RRID: SCR_000432
Seurat	CRAN / Satija Lab	version 4.4.0	Prétraitement, clustering et visualisation du RNA-seq à cellule unique RRID: SCR_016341
shapviz	CRAN	version 0.10.2	Analyse d'interprétabilité du modèle basée sur SHAP