Dywergencja mikrobioty korzeni w różnych siedliskach w oparciu o ważone sieci korelacji

Podsumowanie

Analiza sieci została zastosowana do oceny związku różnych ekologicznych zbiorowisk mikroorganizmów, takich jak gleba, woda i ryzosfera. Przedstawiono protokół dotyczący wykorzystania algorytmu WGCNA do analizy różnych sieci współwystępowania, które mogą wystąpić w zbiorowiskach mikroorganizmów ze względu na różne środowiska ekologiczne.

Streszczenie

Mikrobiom korzeniowy odgrywa ważną rolę we wzroście roślin i adaptacji do środowiska. Analiza sieciowa jest ważnym narzędziem do badania społeczności, które może skutecznie badać relacje interakcji lub model współwystępowania różnych gatunków drobnoustrojów w różnych środowiskach. Celem tego artykułu jest dostarczenie szczegółowych informacji na temat wykorzystania algorytmu ważonej sieci korelacji do analizy różnych sieci współwystępowania, które mogą wystąpić w zbiorowiskach drobnoustrojów ze względu na różne środowiska ekologiczne. Cała analiza eksperymentu jest wykonywana w pakiecie WGCNA. WGCNA to pakiet R do ważonej analizy sieci korelacji. Dane eksperymentalne wykorzystane do zademonstrowania tych metod były danymi dotyczącymi społeczności mikrobiologicznej z bazy danych NCBI (National Center for Biotechnology Information) dla trzech nisz systemu korzeniowego ryżu (Oryza sativa). Wykorzystaliśmy algorytm sieci korelacji ważonej do skonstruowania sieci współliczebności społeczności drobnoustrojów w każdej z trzech nisz. Następnie zidentyfikowano zróżnicowane sieci współobfitości między glebą endosfery, ryzopłaszczyzny i ryzosfery. Ponadto podstawowe rodzaje w sieci uzyskano za pomocą pakietu "WGCNA", który odgrywa ważną rolę regulowaną w funkcjach sieciowych. Metody te umożliwiają naukowcom analizę reakcji sieci mikrobiologicznej na zakłócenia środowiskowe i weryfikację różnych teorii ekologicznej reakcji mikrobiologicznej. Wyniki tych metod wskazują, że znaczące zróżnicowane sieci mikrobiologiczne zidentyfikowane w endosferze, glebie ryzopłaskiej i ryzosferowej ryżu.

Wprowadzenie

Badania mikrobiomu mają ważne implikacje dla zrozumienia i manipulowania procesami ekosystemu^1,2. Populacje drobnoustrojów są połączone wzajemnie oddziałującymi sieciami ekologicznymi, których cechy mogą wpływać na reakcję mikroorganizmów na zmiany środowiskowe^3,4. Co więcej, właściwości tych sieci wpływają na stabilność zbiorowisk mikrobiologicznych i są ściśle związane z funkcją gleby⁵. Analiza sieci ważonej korelacji genów jest obecnie szeroko stosowana w badaniach nad związkiem między genami a społecznościami drobnoustrojów⁶. Wcześniejsze badania koncentrowały się głównie na związkach między sieciami różnych genów lub populacji a światem zewnętrznym⁷. Jednak różnice w sieciach korelacji tworzonych przez populacje drobnoustrojów w różnych warunkach środowiskowych były prawie niezbadane. Celem badań przedstawionych w niniejszym artykule jest dostarczenie spostrzeżeń i szczegółów na temat szybkiej implementacji algorytmu WGCNA do budowy sieci współwystępowania próbek mikrobiomu pobranych w różnych warunkach środowiskowych. Na podstawie wyników analizy oceniliśmy skład i różnice populacji oraz dalej omówiliśmy relacje między różnymi populacjami drobnoustrojów. Zastosowano następujący podstawowy przepływ algorytmu sieci korelacji ważonej⁸. Po pierwsze, należało skonstruować macierz podobieństwa, obliczając współczynnik korelacji Pearsona między profilami ekspresji operacyjnych jednostek taksonomicznych (OTU). Następnie parametry funkcji sąsiedztwa (funkcji przyległości potęgowej lub esicznej) przyjęto za pomocą kryterium topologii bezskalowej, macierz podobieństwa przekształcono w macierz sąsiedztwa, a każda sieć współwystępowania odpowiadała macierzy sąsiedztwa. Użyliśmy hierarchicznego grupowania średnich powiązań w połączeniu z odmiennością opartą na TOM, aby pogrupować OTU ze spójnymi profilami ekspresji w moduły. Następnie obliczyliśmy zależność między konserwatywnymi statystykami a powiązanymi modułami analizy parametrów, ostatecznie identyfikując OTU huba w module. Metody te są szczególnie przydatne do analizy różnic w strukturach sieciowych między różnymi populacjami drobnoustrojów w rozbieżnych warunkach środowiskowych. W niniejszym manuskrypcie szczegółowo opisaliśmy metodę rozwoju sieci koekspresji, analizę różnic między modułami oraz przedstawiliśmy krótki przegląd etapów procedury stosowanej w celu uzyskania gatunków rdzeniowych w różnych sieciach modułowych.

Protokół

1. Pobieranie danych

1. Pobierz dane PRJNA386367 akcesyjnego z bazy NCBI. Z danych PRJNA386367 akcesyjnego wybierz dane dotyczące mikrobiomu ryzosfery, ryzopłaszczyzny i endosfery z roślin ryżu uprawianych przez 14 tygodni na zanurzonym polu ryżowym w Arbuckle w Kalifornii w 2014 r.
   UWAGA: Dane dotyczące mikrobiomu ryzosfery, ryzopłaszczyzny i endosfery zostały przedstawione w tabeli OTU w PRJNA386367 akcesyjnym.

2. Wyznaczanie optymalnej wartości mocy

UWAGA: Pakiet WGCNA zawiera wszystkie następujące parametry funkcjonalne. WGCNA to pakiet R do ważonej analizy sieci korelacji. Kluczowe wiersze poleceń odnoszą się do suplementu S1.

1. W środowisku języka R otwórz oprogramowanie Rstudio i zainstaluj pakiet WGCNA.
2. Załaduj dane i użyj funkcji goodSamplesGenes, aby sprawdzić poprawność danych. Wykonaj wiersze poleceń:
   "gsg = goodSamplesGenes(datExpr0, verbose = 3)
   gsg$allOK "
   Kliknij przycisk Uruchom.
3. Sprawdź, czy nie ma wartości odstających i przechowuj próbki, które spełniają wymagania. Gdy wynik sprawdzenia jest PRAWDZIWY, przejdź do następnego kroku. Zapisz wynik.
4. Funkcja PickSoftThreshold służy do obliczania indeksu bez skalowania R² dwóch grup danych o różnych wartościach mocy. Wykonaj wiersz poleceń:
   "sft = pickSoftThreshold(datExpr0, powerVector = potęgi, pełne = 5)"
   Kliknij przycisk Uruchom.
5. Wizualizacja wyników (Rysunek 1). Wykonaj wiersz poleceń:
   "plot(sft$fitIndices[,1], -sign(sft$fitIndices[,3])*sft$fitIndices[,2],
   xlab="Miękki próg (moc)",ylab="Dopasowanie modelu topologii bez skalowania,podpisane R^2",type="n",
   main = paste("ES_Scale niezależność"));
   text(sft$fitIndeksy[,1], -sign(sft$fitIndeksy[,3])*sft$fitIndeksy[,2],
   labels=moce,cex=cex1,col="czerwony");
   abline(h=0.9,col="czerwony")
   plot(sft$fitIndeksy[,1], sft$fitIndeksy[,5];
   xlab="Miękki próg (moc)",ylab="Średnia łączność", type="n",
   main = paste("ES_Mean łączność"))
   text(sft$fitIndices[,1], sft$fitIndices[,5], labels=powers, cex=cex1,col="red")"
   Kliknij Uruchom.
   UWAGA: Założeniem algorytmu sieci korelacji ważonej jest to, że ustalona struktura sieci współwyrażeń jest zgodna ze standardami kryterium topologii bezskalowej, co zwiększa jej niezawodność. Indeks bez skalowania bliższy 1 wskazuje strukturę sieci, która jest bliższa sieci bez skalowania.
6. Wybierz wartość mocy, gdy wskaźnik bezskalowy R² do kwadratu jest większy niż 0,9 i przejdź do następnego kroku analizy.
   UWAGA: Gdy indeks bez skalowania jest bliski 1, struktura sieci jest bliższa sieci bez skalowania. Analizując dwie lub więcej sieci, konieczne jest podjęcie decyzji, aby każda sieć była zbliżona do wartości mocy sieci bezskalowej, aby zapewnić porównywalność między sieciami o współwyrażanej wartości.

3. Budowa sieci współekspresji i identyfikacja modułu

UWAGA: Na podstawie powyższej obliczonej wartości mocy, budowana jest sieć współwystępowania. Kluczowe wiersze poleceń odnoszą się do suplementu S2.

1. Użyj funkcji sąsiedztwa w pakiecie WGCNA, aby dodać podpisane parametry do budowy symbolicznej sieci współwystępowania. Wykonaj wiersz poleceń:
   "sąsiedztwo = sąsiedztwo(datExpr0, moc = miękkaMoc)"
   Kliknij przycisk Uruchom.
2. Zastosuj funkcję podobieństwa TOM, aby opracować topologiczną sieć nakładającą się i obliczyć sieć odmienności. Wykonaj wiersz poleceń:
   "TOM = TOMpodobieństwo(sąsiedztwo);
   dissTOM = 1-TOM"
   Kliknij Uruchom.
   UWAGA: Parametr signed został dodany w celu ustawienia typu sieci nakładającej się topologii.
3. Użyj funkcji hclust, aby wybrać metodę grupowania hierarchicznego średniego powiązania dla grupowania hierarchicznego. Wykonaj wiersz poleceń:
   "geneTree = hclust(as.dist(dissTOM), metoda = "średnia"); "
   Kliknij przycisk Uruchom.
4. Użyj funkcji cutreeDynamic, aby wykonać dynamiczne wycinanie gałęzi i ustaw parametr minClusterSize na wartość 30. Uzyskaj wynik rozpoznawania modułu. Wykonaj wiersz poleceń:
   "dynamicMods = cutreeDynamic(dendro = geneTree, distM = dissTOM, deepSplit = 2, pamRespectsDendro = FALSE, minClusterSize = minModuleSize); "
   Kliknij Uruchom.
   UWAGA: Minimalny rozmiar modułu nie może być mniejszy niż 30.
5. Oblicz wartość własną każdego modułu OTU za pomocą funkcji modułuEigengenes. Wykonaj wiersz poleceń:
   "MEList = moduleEigengenes(datExpr0, colors = dynamicColors)
   MEs = MEList$eigengenes"
   Kliknij Uruchom.
   UWAGA: Eigen modułu reprezentował ogólny poziom ekspresji OTU w module. Nie była to konkretna OTU, ale pierwszy główny składnik każdego klastra uzyskany w wyniku pojedynczej dekompozycji wartości sieciowej.
6. Wykonaj funkcję klastra na podstawie współczynnika korelacji modułu własnego. Użyj funkcji mergeCloseModules, aby scalić moduły o wartości mniejszej niż 0,25. Wykonaj wiersz poleceń:
   "merge = mergeCloseModules(datExpr0, dynamicColors, cutHeight = MEDissThres, verbose = 3)"
   Kliknij przycisk Uruchom.
7. Na koniec użyj funkcji plotDendroAndColors do wizualizacji, aby uzyskać diagram wyświetlania przypisania modułów dla każdej sieci współwyrażeń (Rysunek 2). Użyj funkcji tabeli, aby wyodrębnić atrybucję modułu odpowiadającą każdemu OT w tabeli przypisania modułu. Wykonaj wiersz poleceń:
   "plotDendroAndColors(geneTree, mergedColors, "Scalony dynamiczny",dendroLabels = FALSE,
   hang = 0,03,addGuide = TRUE, guideHang = 0,05,
   main = "ES_Gene kolory dendrogramu i modułów")"
   Kliknij Uruchom.
   UWAGA: Na diagramie przypisania modułów sieci współwyrażającej różne kolory reprezentują różne moduły, a szary reprezentuje OTU, których nie można zaklasyfikować do żadnego modułu. Większa liczba OTU w szarym module wskazuje, że jakość wstępnego przetwarzania matrycy wyrażeń na wczesnym etapie jest niska.

4. Porównanie modułów

UWAGA: Ta metoda może być użyta do porównania modułów sieciowych dwóch ekologicznych społeczności mikroorganizmów. W tym artykule porównano różnice modułów sieci mikrobiologicznej między endosferą a ryzopłaszczyzną, endosferą a ryzosferą, ryzosferą i ryzopłaszczyzną.

1. Test konserwujący
   1. Załaduj parametry i wyniki dwóch zestawów danych zapisanych w poprzednich krokach.
   2. Ustaw wynik przypisania modułu sieciowego grupy danych mikrobiologicznych jako grupę odniesienia, a drugą grupę jako grupę testową.
   3. Funkcja modulePreservation służy do obliczania wartości parametrów statystycznych konserwatywności Z_summary i medianRank. Wykonaj wiersz poleceń:
      "system.time({mp=modulePreservation(multiExpr,
      multiColor,referenceNetworks=1,
      nPermutation=100, randomSeed=1,quickCor=0,verbose=3)})"
      Kliknij Uruchom.
      UWAGA: Ten wynik może określić ilościowo konserwatywność między modułami. Z_summary>10 oznacza, że dwa moduły są wysoce zachowane, podczas gdy Z_summary<2 oznacza moduły niezachowane. medianRank wyraża względne zachowanie modułu ocenianego przez klasyfikację. Wyższe wartości medianRank oznaczają niezachowane moduły. (Kluczowe wiersze poleceń odnoszą się do Suplementu S3.)
   4. Użyj funkcji wykresu, aby zwizualizować wyniki (Rysunek 3). Pobierz parametry Z_summary i medianRank (tabela 1).
      UWAGA: Moduły sieciowe, które spełniają zarówno wartość Z_summary mniejszą niż 2, jak i medianę wartości rangi na górze, są najbardziej niezachowanymi modułami w dwóch ekologicznych zbiorowiskach mikroorganizmów.
   5. Na podstawie wyników wyżej wymienionych dwóch parametrów statystycznych w celu zidentyfikowania modułu z najbardziej niezachowanym modułem z dwóch sieci.
2. Analiza korelacji członkostwa w module
   1. Ustaw wyniki przypisania modułów dwóch sieci, które zostały ustawione odpowiednio jako odniesienie i grupa testowa.
      UWAGA: Ustawienia muszą być takie same jak w teście konserwacji.
   2. Użyj funkcji corPvalueStudent, aby wyodrębnić wartość kME (członkostwo w module) dla każdego OTU w kilku modułach kandydujących.
      Wykonaj wiersz poleceń:
      "Pvalue = as.data.frame(corPvalueStudent(as.matrix
      (ModuleMembership), próbki))"
      Kliknij Uruchom.
      UWAGA: kME oznacza stopień członkostwa w module. ME oznacza moduł eigen, który reprezentuje ogólny poziom ekspresji OTU w module. kME jest współczynnikiem korelacji między każdą OTU a ME. Określ ilościowo znaczenie OTU w sieci na podstawie wartości kME OTU. (Kluczowe wiersze poleceń odnoszą się do Suplementu S4.)
   3. Następnie użyj funkcji verboseScatterplot, aby obliczyć współczynnik korelacji wartości kME odpowiednich OTU w dwóch sieciach i narysuj diagram analizy korelacji (Rysunek 4).
      Wykonaj wiersz poleceń:
      "verboseScatterplot(abs(TModuleMembership
      [TmoduleGenes, Tcolumn]),
      abs(NModuleMembership[NmoduleGenes, Ncolumn]),
      xlab = paste("kME in", "ES"),
      ylab = paste("kME in", "RP"),
      main = paste("jasnożółty"),
      cex.main = 1.7, cex.lab = 1.6, cex.axis = 1.6, col = modulecolor)"
      Kliknij przycisk Uruchom.
   4. Wybierz moduł o najmniejszym współczynniku korelacji wartości kME OTU obu sieci. Weź pod uwagę, że ten moduł ma największą różnicę z tych dwóch sieci.

5. Analiza modułu mikrobiologicznej sieci różnicowej

1. Uzyskaj dane dotyczące dominującej gromady bakterii poprzez analizę statystyczną zestawu sekwencji OTU modułu o największej różnicy.
   UWAGA: Zestaw sekwencji OTU modułu z największą różnicą jest sumowany przez taksonomię gromady. Dominująca gromada bakterii stanowiła ponad 10%.
2. Następnie należy użyć funkcji exportNetworkToCytoscape, aby uzyskać plik zawierający informacje o relacji interakcji OTU w największym module różnicowym.
   Wykonaj wiersz poleceń:
   "cyt = exportNetworkToCytoscape(modTOM,
   edgeFile = paste("NOWY-ES_CytoscapeInput-krawędzie-", moduły , ".txt", sep=""),
   nodeFile = paste("NOWY-ES_CytoscapeInput-węzły-", moduły, ".txt", sep=""),
   weighted = TRUE, threshold = 0.5, nodeNames = modProbes,
   altNodeNames = modGenes, nodeAttr = moduleColors[wModule])"
   Kliknij przycisk Uruchom.
3. Zaimportuj plik do Cytoscape. Ustaw próg na 0,5 i dostosuj inne parametry zgodnie z potrzebami.
4. Skonstruuj sieć współwystępowania mikroorganizmów różnicowych (Rysunek 5).
5. Uzyskano informacje o podstawowym rodzaju, który pełni najważniejszą rolę regulacyjną w sieci.
   UWAGA: Zgodnie z wartością kME OUT można zdefiniować rdzeń rodzaju.
6. Na koniec oceniono funkcje rodzaju rdzeniowego i przeanalizowano jego wpływ na całą sieć różnic.

Wyniki

Reprezentatywne wyniki w tym artykule zostały pobrane z danych mikrobiomu korzenia ryżu California Abaker z 2014 roku w bazie danych NCBI (PRJNA386367)⁹. Dane obejmują próbki mikrobiomu ryzosfery, ryzopłaszczyzny i endosfery z roślin ryżu uprawianych przez 14 tygodni na zanurzonym polu ryżowym. Użyliśmy algorytmu WGCNA, aby wybrać wartość mocy, która spełniała trzy sieci, które były zbliżone do sieci bezskalowej (Rysunek 1) i opracowaliśmy trzy sieci koekspresji (

Dyskusja

Sieci korelacyjne są coraz częściej wykorzystywane w zastosowaniach bioinformatycznych. WGCNA to metoda biologii systemowej służąca do opisowej analizy relacji między różnymi elementami systemu biologicznego¹². Pakiet oprogramowania R był używany we wcześniejszych pracach nad WGCNA 13,14,15. Pakiet obejmuje funkcje do budowy sieci, wykrywania modułów, obliczeń właściwości topologicznych, symulacji danych,...

Materiały

Lista materiałów użytych w tym artykule

Nazwa	Firma	Numer katalogowy	Komentarze
R	University of Auckland	w wersji 4.0.2	R to wolne środowisko oprogramowania do obliczeń statystycznych i grafiki. Kompiluje się i działa na wielu różnych platformach UNIX, Windows i MacOS.
RStdio	JJ Allaire	w wersji 1.4.1103	RStudio IDE to zestaw zintegrowanych narzędzi zaprojektowanych w celu zwiększenia produktywności w językach R i Python.
Cytoscape		w wersji 3.7.1	Cytoscape to platforma oprogramowania typu open source do wizualizacji złożonych sieci i integrowania ich z dowolnym typem danych atrybutów.
Baza danych			Narodowe Centrum Informacji Biotechnologicznej przyczynia się do rozwoju nauki i zdrowia, zapewniając dostęp do informacji biomedycznych i genomicznych.