Wstępna analiza i walidacja danych sekwencjonowania CUT&RUN

Podsumowanie

Ten protokół prowadzi początkujących bioinformatyków przez wprowadzający potok analizy CUT&RUN, który umożliwia użytkownikom przeprowadzenie wstępnej analizy i walidacji danych sekwencjonowania CUT&RUN. Wykonanie opisanych tutaj kroków analizy, w połączeniu z adnotacją pików, pozwoli użytkownikom na uzyskanie mechanistycznego wglądu w regulację chromatyny.

Streszczenie

Technika CUT&RUN ułatwia wykrywanie interakcji białko-DNA w całym genomie. Typowe zastosowania CUT&RUN obejmują profilowanie zmian w modyfikacjach ogona histonowego lub mapowanie zajętości chromatyny czynnika transkrypcyjnego. Powszechne zastosowanie CUT&RUN jest częściowo napędzane przez zalety techniczne w porównaniu z konwencjonalnym sekwencją ChIP, które obejmują niższe wymagania dotyczące wprowadzania komórek, niższe wymagania dotyczące głębokości sekwencjonowania oraz zwiększoną czułość przy zmniejszonym sygnale tła ze względu na brak środków sieciujących, które w przeciwnym razie maskują epitopy przeciwciał. Powszechne przyjęcie CUT&RUN zostało również osiągnięte dzięki hojnemu udostępnianiu odczynników przez laboratorium Henikoff oraz opracowaniu komercyjnych zestawów w celu przyspieszenia wdrożenia przez początkujących. Wraz ze wzrostem technicznego zastosowania CUT&RUN, analiza sekwencjonowania CUT&RUN i walidacja stają się krytycznymi wąskimi gardłami, które muszą zostać pokonane, aby umożliwić pełne przyjęcie przez zespoły laboratoryjne, w których przeważają mokre laboratoria. Analiza CUT&RUN zazwyczaj rozpoczyna się od kontroli jakości surowych odczytów sekwencjonowania w celu oceny głębokości sekwencjonowania, jakości odczytu i potencjalnych odchyleń. Odczyty są następnie dopasowywane do referencyjnego zestawu sekwencji genomu, a następnie wykorzystuje się kilka narzędzi bioinformatycznych do opisywania regionów genomu wzbogacenia białek, potwierdzania możliwości interpretacji danych i wyciągania wniosków biologicznych. Chociaż opracowano wiele potoków analizy in silico w celu wsparcia analizy danych CUT&RUN, ich złożona, wielomodułowa struktura i wykorzystanie wielu języków programowania sprawiają, że platformy te są trudne dla początkujących bioinformatyków, którzy mogą nie znać wielu języków programowania, ale chcą zrozumieć procedurę analizy CUT&RUN i dostosować swoje potoki analityczne. W tym miejscu udostępniamy jednojęzyczny protokół potoku analizy krok po kroku CUT&RUN, przeznaczony dla użytkowników o dowolnym poziomie doświadczenia bioinformatycznego. Protokół ten obejmuje przeprowadzanie krytycznych kontroli jakości w celu sprawdzenia, czy dane sekwencjonowania są odpowiednie do interpretacji biologicznej. Oczekujemy, że przestrzeganie protokołu wprowadzającego przedstawionego w tym artykule w połączeniu z adnotacją o pikach w dół pozwoli użytkownikom na wyciąganie wniosków biologicznych z własnych zestawów danych CUT&RUN.

Wprowadzenie

Możliwość mierzenia interakcji między białkami a genomowym DNA jest fundamentalna dla zrozumienia biologii regulacji chromatyny. Skuteczne testy, które mierzą zajętość chromatyny dla danego białka, dostarczają co najmniej dwóch kluczowych informacji: i) lokalizacji genomowej oraz ii) obfitości białka w danym regionie genomu. Śledzenie zmian rekrutacji i lokalizacji białka będącego przedmiotem zainteresowania w chromatynie może ujawnić bezpośrednie docelowe loci białka i ujawnić mechanistyczne role tego białka w procesach biologicznych opartych na chromatynie, takich jak regulacja transkrypcji, naprawa DNA lub replikacja DNA. Dostępne obecnie techniki profilowania interakcji białko-DNA umożliwiają naukowcom badanie regulacji w niespotykanej dotąd rozdzielczości. Taki postęp techniczny był możliwy dzięki wprowadzeniu nowych technik profilowania chromatyny, które obejmują opracowanie przez laboratorium Henikoff metody rozszczepienia pod celami i uwalniania za pomocą nukleazy (Cleavage Under Targets and Release Using Nucleaase). CUT&RUN oferuje kilka zalet technicznych w porównaniu z konwencjonalną immunoprecypitacją chromatyny (ChIP), które obejmują niższe wymagania dotyczące wprowadzania komórek, niższe wymagania dotyczące głębokości sekwencjonowania oraz zwiększoną czułość przy zmniejszonym sygnale tła ze względu na brak środków sieciujących, które w przeciwnym razie maskują epitopy przeciwciał. Przyjęcie tej techniki do badania regulacji chromatyny wymaga dogłębnego zrozumienia zasady leżącej u podstaw tej techniki oraz zrozumienia, jak analizować, walidować i interpretować dane CUT&RUN.

Procedura CUT&RUN rozpoczyna się od wiązania komórek z konkanawaliną A sprzężoną z kulkami magnetycznymi, aby umożliwić manipulowanie niską liczbą komórek podczas całej procedury. Wyizolowane komórki są przepuszczane przy użyciu łagodnego detergentu, aby ułatwić wprowadzenie przeciwciała, które jest skierowane przeciwko białemu zainteresowaniu. Nukleaza mikrokokkowa (MNaza) jest następnie rekrutowana do związanego przeciwciała za pomocą znacznika białka A lub białka A / G połączonego z enzymem. Wapń jest wprowadzany w celu zainicjowania aktywności enzymatycznej. W wyniku trawienia MNazy powstają mononukleosomalne kompleksy DNA-białko. Wapń jest następnie chelatowany w celu zakończenia reakcji trawienia, a krótkie fragmenty DNA z trawienia MNazy są uwalniane z jąder, a następnie poddawane oczyszczaniu DNA, przygotowaniu biblioteki i sekwencjonowaniu o wysokiej przepustowości¹ (Rysunek 1).

Podejścia in silico do mapowania i ilościowego określania obecności białek w całym genomie rozwinęły się równolegle z metodami mokrego laboratorium, używanymi do wzbogacania tych interakcji DNA-białko. Identyfikacja regionów wzbogaconych sygnałów (pików) jest jednym z najbardziej krytycznych etapów analizy bioinformatycznej. Początkowe metody analizy ChIP-seq wykorzystywały algorytmy, takie jak MACS² i SICER³, które wykorzystywały modele statystyczne do odróżnienia bona fide miejsc wiązania białko-DNA od szumu tła. Jednak niższy szum tła i wyższa rozdzielczość danych CUT&RUN sprawiają, że niektóre programy wywołujące wartości szczytowe stosowane w analizie ChIP-seq nie nadają się do analizy CUT&RUN⁴. Wyzwanie to uwypukla potrzebę opracowania nowych narzędzi, które lepiej nadają się do analizy danych CUT&RUN. SEACR⁴ reprezentuje jedno z takich narzędzi, które zostało ostatnio opracowane w celu umożliwienia wywoływania szczytów z danych CUT&RUN, przy jednoczesnym pokonywaniu ograniczeń związanych z narzędziami zwykle stosowanymi do analizy ChIP-seq.

Biologiczne interpretacje z danych sekwencjonowania CUT&RUN są pobierane z danych wyjściowych poniżej wywołania piku w potoku analizy. Można zaimplementować kilka funkcjonalnych programów adnotacyjnych w celu przewidywania potencjalnego znaczenia biologicznego wywoływanych pików na podstawie danych CUT&RUN. Na przykład projekt Gene Ontology (GO) zapewnia ugruntowaną funkcjonalną identyfikację genów będących przedmiotem zainteresowania^5,6,7. Różne narzędzia i zasoby programowe ułatwiają analizę GO w celu ujawnienia genów i zestawów genów wzbogaconych wśród szczytów CUT&RUN^{8,9,10,11,12,13,14}. Ponadto oprogramowanie do wizualizacji, takie jak Deeptools¹⁵, Integrative genomics viewer (IGV)¹⁶ i UCSC Genome Browser¹⁷ umożliwiają wizualizację rozkładu sygnałów i wzorców w interesujących regionach genomu.

Zdolność do wyciągania interpretacji biologicznych z danych CUT&RUN zależy przede wszystkim od walidacji jakości danych. Kluczowe komponenty do walidacji obejmują ocenę: i) jakości sekwencjonowania biblioteki CUT&RUN, ii) podobieństwa replikacji oraz iii) dystrybucji sygnału w centrach szczytów. Zakończenie walidacji wszystkich trzech komponentów ma kluczowe znaczenie dla zapewnienia wiarygodności próbek z biblioteki CUT&RUN i wyników dalszej analizy. W związku z tym konieczne jest opracowanie wstępnych przewodników po analizie CUT&RUN, aby umożliwić początkującym bioinformatykom i badaczom w laboratoriach mokrych przeprowadzenie takich etapów walidacji w ramach standardowych procesów analizy CUT&RUN.

Wraz z rozwojem eksperymentu CUT&RUN w laboratorium mokrym, różne metody analizy in silico CUT&RUN, takie jak CUT&RUNTools 2.0^18,19, nf-core/cutandrun²⁰ oraz CnRAP²¹, zostały opracowane do obsługi analizy danych CUT&RUN. Narzędzia te zapewniają zaawansowane podejścia do analizy jednokomórkowych i zbiorczych zestawów danych CUT&RUN i CUT&Tag. Jednak stosunkowo złożona modułowa struktura programu i wymagana znajomość wielu języków programowania do przeprowadzania tych potoków analitycznych może utrudniać przyjęcie ich przez początkujących bioinformatyków, którzy chcą dokładnie zrozumieć etapy analizy CUT&RUN i dostosować własne potoki. Obejście tej bariery wymaga nowego, wprowadzającego potoku analizy CUT&RUN, który jest dostarczany w postaci prostych skryptów krok po kroku zakodowanych przy użyciu prostego, pojedynczego języka programowania.

W tym artykule opisujemy prosty, jednojęzyczny protokół potoku analizy CUT&RUN, który dostarcza krok po kroku skrypty wspierane szczegółowymi opisami, aby umożliwić nowym i początkującym użytkownikom przeprowadzanie analizy sekwencjonowania CUT&RUN. Programy używane w tym potoku są publicznie dostępne przez oryginalne grupy deweloperów. Główne kroki opisane w tym protokole obejmują wyrównanie odczytu, wywołanie pików, analizę funkcjonalną i, co najważniejsze, etapy walidacji w celu oceny jakości próbki w celu określenia przydatności i wiarygodności danych do interpretacji biologicznej (Rysunek 2). Co więcej, potok ten daje użytkownikom możliwość porównywania wyników analiz z publicznie dostępnymi zestawami danych CUT&RUN. Ostatecznie, ten protokół procesu analizy CUT&RUN służy jako przewodnik wprowadzający i odniesienie dla początkujących analityków bioinformatycznych i badaczy w laboratoriach mokrych.

Protokół

UWAGA: Informacje dotyczące plików fastq dla metody CUT&RUN w zbiorze GSE126612 znajdują się w Tabeli 1. Informacje dotyczące aplikacji programistycznych użytych w niniejszym badaniu wymieniono w Tabeli Materiałów.

1. Pobieranie potoku Easy-Shells_CUTnRUN z jego strony na Githubie

1. Otwórz terminal w systemie operacyjnym.
   UWAGA: Jeśli użytkownik nie wie, jak otworzyć terminal w systemach macOS i Windows, należy zapoznać się z tą stroną (https://discovery.cs.illinois.edu/guides/System-Setup/terminal/). W przypadku systemu Linux należy zapoznać się z tą stroną (https://www.geeksforgeeks.org/how-to-open-terminal-in-linux/).
2. Pobierz skompresowany potok analizy z serwisu Github, wpisując w terminalu wget https://github.com/JunwooLee89/Easy-Shells_CUTnRUN/archive/refs/heads/main.zip -O ~/Desktop/Easy-Shells_CUTnRUN.zip .
3. Po pobraniu pliku zip, rozpakuj go, wpisując w terminalu unzip ~/Desktop/Easy-Shells_CUTnRUN.zip -d ~/Desktop/ .
4. Po rozpakowaniu usuń plik zip, wpisując w terminalu rm ~/Desktop/Easy-Shells_CUTnRUN.zip i zmień nazwę folderu, wpisując mv ~/Desktop/Easy-Shells_CUTnRUN-master ~/Desktop/Easy-Shells_CUTnRUN.
5. Po usunięciu skompresowanego pliku, wpisz w terminalu chmod +x ~/Desktop/Easy-Shells_CUTnRUN/script/*.sh, aby nadać uprawnienia do wykonywania wszystkich skryptów powłoki w katalogu roboczym. Od teraz wystarczy wpisać ścieżkę i nazwę tych skryptów w terminalu lub przeciągnąć skrypty do terminala i nacisnąć enter, aby je uruchomić.
   UWAGA: Powłoka Bash jest zazwyczaj wstępnie zainstalowana w większości dystrybucji Linuxa. Jednak nowsze wersje macOS nie oferują już wstępnie zainstalowanej powłoki Bash. Jeśli system nie posiada powłoki Bash, należy ją najpierw zainstalować. Instrukcje dotyczące instalacji powłoki Bash w systemie Linux (https://ioflood.com/blog/install-bash-shell-linux/) oraz macOS (https://www.cs.cornell.edu/courses/cs2043/2024sp/styled-3/#:~:text=The%20first%20thing%20you%20will,you%20will%20see%20the%20following:) znajdują się pod poniższymi linkami. Te krokowe skrypty powłoki zostały napisane tak, aby utworzyć jeden folder ~/Desktop/GSE126612 i przeprowadzić większość analizy CUT&RUN w tym katalogu bez potrzeby wprowadzania jakichkolwiek modyfikacji. Jeśli użytkownik rozumie sposób działania tych skryptów powłoki, może je zrewidować i dostosować do analizy innych zbiorów danych CUT&RUN oraz zmodyfikować opcje zgodnie ze specyficznymi potrzebami projektu. Do odczytu i edycji tych skryptów powłoki można wykorzystać Visual Studio Code (https://code.visualstudio.com/) jako jedną z opcji łatwego w użyciu programu dostępnego dla głównych systemów operacyjnych.

2. Instalacja programów wymaganych dla Easy Shells CUTnRUN

1. Spośród skryptów powłoki o nazwie Script_01_installation_***.sh, znajdź skrypt, którego nazwa zawiera typ systemu operacyjnego używanego przez użytkownika. Obecnie Easy Shells CUTnRUN obsługuje skrypty instalacyjne dla systemów macOS, Debian/Ubuntu oraz systemów opartych na CentOS/RPM.
2. Otwórz terminal i wpisz echo $SHELL, aby sprawdzić domyślną powłokę w aktywnym terminalu. Jeśli domyślną powłoką w bieżącym terminalu jest Bash, użytkownicy mogą zobaczyć następujący komunikat: /path/to/bash (lub podobny, np. /bin/bash) w terminalu.
3. Jeśli domyślną powłoką nie jest Bash, ustaw powłokę Bash jako domyślną, wpisując w terminalu chsh -s $(which bash) . Jeśli terminal używa powłoki Bash jako domyślnej, pomiń ten krok.
4. W terminalu uruchom instalacyjny skrypt powłoki, wpisując ~/Desktop/Easy-Shells_CUTnRUN/scripts/Script_01_installation_***.sh lub przeciągając plik skryptu do terminala i zatwierdzając klawiszem enter.
5. Przeczytaj plik Test_README.md w folderze /path/to/SEACR-1.3/Testfiles. Postępuj zgodnie z instrukcjami zawartymi w pliku README, aby upewnić się, że program SEACR w systemie użytkownika działa prawidłowo.
   UWAGA: W celu uzyskania prawidłowych wyników wywoływania pików (peak calling) z danych CUT&RUN kluczowe jest zwalidowanie funkcji SEACR za pomocą plików testowych udostępnionych na stronie SEACR na Githubie. Dlatego niezwłocznie po instalacji SEACR postępuj zgodnie z instrukcjami zawartymi w  Test_README.md w /path/to/SEACR-1.3/Testfiles. Mimo że Easy Shells CUTnRUN dostarcza skrypty instalacyjne dla niektórych systemów operacyjnych, skrypty te mogą nie zadziałać w systemach niektórych użytkowników podczas instalacji wszystkich programów wymaganych dla Easy Shells CUTnRUN. W przypadku wystąpienia problemów z instalacją zapoznaj się z oryginalną stroną programu, którego nie udało się zainstalować, lub poproś o pomoc za pomocą strony z problemami (issues) Easy Shells CUTnRUN na githubie (https://github.com/JunwooLee89/Easy-Shells_CUTnRUN/issues).

3. Pobieranie publicznie dostępnego zbioru danych CUT&RUN z Sequence Read Archive (SRA)

1. Otwórz terminal i wpisz echo $SHELL, aby sprawdzić domyślną powłokę w aktywnym terminalu. Jeśli powłoką domyślną w aktualnym terminalu jest Bash, użytkownicy zobaczą w terminalu następującą ścieżkę: /path/to/bash (lub podobny komunikat, np. /bin/bash).
2. Jeśli powłoka domyślna nie jest powłoką Bash, ustaw Bash jako powłokę domyślną, wpisując w terminalu chsh -s $(which bash) . Jeśli terminal używa powłoki Bash jako domyślnej, pomiń ten krok.
3. Wpisz ~/Desktop/Easy-Shells_CUTnRUN/scripts/Script_02_download-fastq.sh w terminalu lub przeciągnij plik skryptu powłoki do terminala i naciśnij enter.
   UWAGA: Ten skrypt spowoduje: (i) Utworzenie jednego folderu (~/Desktop/GSE126612/fastq) i pobranie listy plików SRA zapisanej w pliku tekstowym (~/Desktop/Easy-Shells_CUTnRUN/sample_info/SRR_list.txt) do folderu fastq. Przykładowo, SRR_list.txt zawiera pliki fastq dla podzbioru próbek CUT&RUN z serii GSE126612. (ii) Pobranie surowych plików fastq do folderu fastq. (iii) Utworzenie jednego folderu (~/Desktop/GSE126612/log/fastq) i zapisanie w nim pliku dziennika (download-fastq_log.txt) oraz pliku z informacjami o pobranych próbkach (SRR_list_info.txt).
4. Po uruchomieniu skryptu sprawdź plik dziennika. Jeśli w pliku dziennika pojawi się jakikolwiek komunikat o błędzie, napraw błąd i ponownie wykonaj krok 3.3. W przypadku problemów z rozwiązaniem błędu, poproś o pomoc na stronie problemów GitHub Easy Shells CUTnRUN (https://github.com/JunwooLee89/Easy-Shells_CUTnRUN/issues).
   UWAGA: Aby ułatwić ćwiczenie tego rurociągu analizy CUT&RUN, z SRA pobrano następujące ogólnodostępne próbki: jedną próbkę z kontroli pozorowanej (IgG), trzy próbki białka architektury chromatyny i czynnika transkrypcyjnego (CTCF), cztery próbki odpowiadające „aktywnej” modyfikacji histonowej (H3K27Ac) oraz trzy próbki odpowiadające regionom inicjacji transkrypcji znakowanym przez polimerazę RNA II (RNAPII-S5P). Sekwencjonowanie wykonano w trybie paired-end, dlatego dla każdej próbki sparowane są dwa pliki.

4. Wstępna kontrola jakości surowych plików z sekwencjonowania

1. Otwórz terminal i wpisz echo $SHELL , aby sprawdzić domyślną powłokę w aktywnym terminalu. Jeśli domyślną powłoką w bieżącym terminalu jest Bash, użytkownicy zobaczą w terminalu następujący komunikat: /path/to/bash (lub podobny, np. /bin/bash).
2. Jeśli domyślna powłoka nie jest Bashem, ustaw powłokę Bash jako domyślną, wpisując w terminalu chsh -s $(which bash). Jeśli terminal używa powłoki Bash jako domyślnej, pomiń ten krok.
3. Wpisz w terminalu ~/Desktop/Easy-Shells_CUTnRUN/scripts/Script_03_fastQC.sh lub przeciągnij skrypt powłoki do terminala i zatwierdź klawiszem enter.
   UWAGA: Ten skrypt powłoki: (i) uruchomi program FastQC dla wszystkich surowych plików fastq w folderze ~/Desktop/GSE126612/fastq i zapisze raporty kontroli jakości w folderze ~/Desktop/GSE126612/fastqc.1st . (ii) utworzy plik dziennika (fastqc.1st.log.SRR-number.txt) dla każdego uruchomienia FastQC w folderze z logami (~/Desktop/GSE126612/log/fastqc.1st).
4. Po zakończeniu działania skryptu powłoki przejrzyj plik dziennika, aby potwierdzić pomyślne wykonanie procesu. Jeśli w pliku dziennika pojawi się jakikolwiek komunikat o błędzie, napraw błąd i powtórz krok 4.3. W przypadku problemów z rozwiązaniem błędu poproś o pomoc za pomocą strony problemów (issues) projektu Easy Shells CUTnRUN w serwisie GitHub (https://github.com/JunwooLee89/Easy-Shells_CUTnRUN/issues).
   UWAGA: Wśród plików wynikowych pliki fastqc.html zawierają przejrzyste wyniki kontroli jakości. W przypadku wystąpienia poważnych problemów z jakością skonsultuj się z bioinformatykami, aby określić przydatność danych do dalszych analiz. Podobne raporty kontroli jakości są wykorzystywane do potwierdzenia poprawy jakości danych po przycinaniu adapterów. Aby użyć tego skryptu dla innych zestawów danych, edytuj ścieżki katalogów roboczych i wyjściowych zgodnie z potrzebami użytkownika. Istotną różnicą przy interpretacji QC odczytów CUT&RUN w porównaniu do ChIP-seq jest to, że zduplikowane odczyty w CUT&RUN niekoniecznie wskazują na duplikaty PCR. Wynika to z faktu, że zrekrutowana nukleaza MNase trawi DNA w tych samych lub podobnych miejscach w grupach eksperymentalnych.

5. Kontrola jakości i przycinanie adapterów dla surowych plików z sekwencjonowania

1. Otwórz terminal i wpisz echo $SHELL, aby sprawdzić domyślną powłokę w aktywnym terminalu. Jeśli domyślną powłoką w bieżącym terminalu jest Bash, użytkownicy mogą zobaczyć następujący komunikat: /path/to/bash (lub podobną informację, np. /bin/bash) w terminalu.
2. Jeśli domyślną powłoką nie jest Bash, ustaw powłokę Bash jako domyślną, wpisując w terminalu chsh -s $(which bash). Jeśli terminal używa powłoki Bash jako domyślnej, pomiń ten krok.
3. Wpisz ~/Desktop/Easy-Shells_CUTnRUN/scripts/Script_04_trimming.sh w terminalu lub przeciągnij skrypt Script_04_trimming.sh do terminala i zatwierdź klawiszem Enter.
   UWAGA: Ten skrypt powłoki będzie: (i) uruchamiać program Trim-Galore dla wszystkich surowych plików fastq w ~/Desktop/GSE126612/fastq, aby przeprowadzić przycinanie adapterów i kontrolę jakości. (ii) tworzyć jeden folder (~/Desktop/GSE126612/trimmed) i zapisywać w nim pliki wyjściowe programu Trim-Galore. (iii) tworzyć jeden folder na logi (~/Desktop/GSE126612/log/trim_galore) i zapisywać plik logu trim_galore_log_RSS-number.txt dla każdego uruchomienia programu Trim-Galore.
4. Po zakończeniu procesu dokładnie przejrzyj plik logu. Jeśli w pliku logu pojawi się jakikolwiek komunikat o błędzie, napraw błąd i powtórz krok 5.3. W przypadku problemów z rozwiązaniem kwestii, poproś o pomoc za pośrednictwem strony z problemami (issues) projektu Easy Shells CUTnRUN na GitHubie (https://github.com/JunwooLee89/Easy-Shells_CUTnRUN/issues). 
5. Po zakończeniu tego procesu porównaj wyjściowe pliki .html z plikami fastqc.html utworzonymi w kroku 4.3. Zmień ścieżki katalogów wejściowych i wyjściowych, aby przeprowadzić krok przycinania dla plików fastq znajdujących się w innych lokalizacjach.

6. Pobieranie indeksu bowtie2 dla genomów referencyjnych dla próbek rzeczywistych i kontroli spike-in

1. Otwórz terminal i wpisz echo $SHELL , aby sprawdzić domyślną powłokę w aktywnym terminalu. Jeśli domyślną powłoką w bieżącym terminalu jest Bash, użytkownicy mogą zobaczyć w terminalu następujący komunikat: /path/to/bash (lub podobny, np. /bin/bash).
2. Jeśli domyślna powłoka nie jest powłoką Bash, ustaw Bash jako domyślną powłokę, wpisując w terminalu chsh -s $(which bash). Jeśli terminal używa powłoki Bash jako domyślnej, pomiń ten krok.
3. Wpisz w terminalu ~/Desktop/Easy-Shells_CUTnRUN/scripts/Script_05_bowtie2-index.sh lub przeciągnij skrypt powłoki do terminala i zatwierdź klawiszem enter.
   UWAGA: Ten skrypt spowoduje: (i) pobranie indeksów Bowtie2 dla rzeczywistych genomów referencyjnych próbek (człowiek; hg19; użyty w oryginalnej publikacji²²) oraz genomów referencyjnych kontroli Spike-in (drożdże pączkujące; R64-1-1) do folderu bowtie2-index (~/Desktop/Easy-Shells_CUTnRUN/bowtie2-index). (iii) zapisanie pliku dziennika (bowtie2-index-log.txt) w katalogu logów (~/Desktop/GSE126612/log/bowtie2-index).
4. Po zakończeniu uruchamiania sprawdź plik dziennika. W przypadku wystąpienia jakiegokolwiek komunikatu o błędzie, napraw błąd i powtórz krok 6.3. Jeśli pojawi się problem z rozwiązaniem kwestii, poproś o pomoc za pośrednictwem strony problemów GitHub dla Easy Shells CUTnRUN (https://github.com/JunwooLee89/Easy-Shells_CUTnRUN/issues). 
   UWAGA: Obecnie indeksy Bowtie2 dla różnych genomów referencyjnych są dostępne na stronie Bowtie2 (https://bowtie-bio.sourceforge.net/bowtie2/manual.shtml). Użytkownicy mogą edytować Script_05_bowtie2-index.sh, aby pobrać dowolny indeks Bowtie2 zgodnie z ich wymaganiami. Jeśli użytkownik nie może zlokalizować indeksu Bowtie2 dla interesującego go genomu referencyjnego, należy odnaleźć pliki fasta sekwencji genomu referencyjnego w:
   1. FTP Ensembl (https://ftp.ensembl.org/pub/current_fasta/)
   2. na stronie internetowej UCSC (https://hgdownload.soe.ucsc.edu/downloads.html) 
   3. lub w innych bazach danych specyficznych dla danego gatunku.
      Po zlokalizowaniu plików fasta sekwencji genomu referencyjnego, utwórz indeks Bowtie2 dla pobranego genomu referencyjnego, postępując zgodnie z sekcją "The bowtie2-build indexer" (https://bowtie-bio.sourceforge.net/bowtie2/manual.shtml#the-bowtie2-build-indexer) na stronie internetowej Bowtie2.

7. Mapowanie przyciętych odczytów sekwencjonowania CUT&RUN do genomów referencyjnych

1. Otwórz terminal i wpisz echo $SHELL , aby sprawdzić domyślną powłokę w aktywnym terminalu. Jeśli domyślną powłoką w aktualnym terminalu jest Bash, użytkownicy zobaczą w terminalu następujący komunikat: /path/to/bash (lub podobny, np. /bin/bash).
2. Jeśli domyślna powłoka to nie Bash, ustaw powłokę Bash jako domyślną, wpisując w terminalu chsh -s $(which bash). Jeśli terminal używa powłoki Bash jako domyślnej, pomiń ten krok.
3. Wpisz w terminalu ~/Desktop/Easy-Shells_CUTnRUN/scripts/Script_06_bowtie2-mapping.sh lub przeciągnij plik skryptu powłoki do terminala i zatwierdź klawiszem enter.
   UWAGA: Ten skrypt powłoki wykona następujące czynności: (i) uruchomi program bowtie2 w celu niezależnego zmapowania wszystkich przyciętych pod kątem adapterów i jakości plików fastq na genomy referencyjne zarówno eksperymentalne (ludzki; hg19), jak i kontrolne spike-in (drożdże pączkowe; R64-1-1). (ii) uruchomi funkcję samtools view, aby skompresować pliki zmapowanych par odczytów do formatu bam. (iii) utworzy folder (~/Desktop/GSE126612/bowtie2-mapped) i zapisze w nim skompresowane pliki zmapowanych par odczytów. (iv) utworzy folder (~/Desktop/GSE126612/log/bowtie2-mapped) i zapisze logi procesu mapowania w formie plików tekstowych bowtie2_log_hg19_SRR-number.txt dla par odczytów zmapowanych na genom referencyjny hg19 oraz bowtie2_log_R64-1-1_SRR-number.txt dla par odczytów zmapowanych na R64-1-1, aby wskazać wydajność mapowania w folderze logów bowtie2-mapping.
4. Po zakończeniu uruchamiania sprawdź plik logu. Jeśli w pliku logu pojawi się jakikolwiek komunikat o błędzie, popraw błąd i uruchom skrypt powłoki ponownie. W przypadku problemów z rozwiązaniem kwestii, poproś o pomoc, korzystając ze strony zgłoszeń GitHub Easy Shells CUTnRUN (https://github.com/JunwooLee89/Easy-Shells_CUTnRUN/issues).
   UWAGA: Ten skrypt powłoki uruchamia bowtie2 z opcjami mapowania plików sekwencjonowania paired-end w celu znalezienia zgodnych zmapowanych par odczytów o długości fragmentów od 10 bp do 700 bp. Opisy opcji można znaleźć, wpisując bowtie2 --help w terminalu lub odwiedzając stronę internetową bowtie2 (https://bowtie-bio.sourceforge.net/bowtie2/manual.shtml#the-bowtie2-aligner), aby zrozumieć i w razie potrzeby zmienić opcje. Użyj tego skryptu do mapowania innych plików fastq, zmieniając ścieżkę i format nazwy plików fastq oraz indeksów Bowtie2.

8. Sortowanie i filtrowanie plików z zmapowanymi parami odczytów

1. Otwórz terminal i wpisz echo $SHELL , aby sprawdzić domyślną powłokę w aktywnym terminalu. Jeśli domyślną powłoką w obecnym terminalu jest Bash, użytkownicy zobaczą w terminalu następującą treść: /path/to/bash (lub podobny komunikat, np. /bin/bash).
2. Jeśli domyślna powłoka nie jest Bashem, ustaw powłokę Bash jako domyślną, wpisując w terminalu „chsh -s $(which bash)”. Jeśli terminal używa powłoki Bash jako domyślnej, pomiń ten krok.
3. Wpisz w terminalu ~/Desktop/Easy-Shells_CUTnRUN/scripts/Script_07_filter-sort-bam.sh lub przeciągnij plik skryptu powłoki do terminala i naciśnij enter.
   UWAGA: Ten skrypt wykona następujące czynności: (i) Uruchomi funkcję samtools view dla wszystkich skompresowanych plików z mapowanymi parami odczytów w folderze ~/Desktop/GSE126612/bowtie2-mapped, aby odfiltrować pary odczytów zmapowane do niekanonicznych regionów chromosomowych, publicznie adnotowanych regionów czarnej listy (blacklist) oraz regionów powtórzeń TA. (ii) Wykona funkcję samtools sort, aby posortować odfiltrowane pliki bam według nazw fragmentów lub współrzędnych w tym samym katalogu. (iii) Utworzy plik dziennika (log) dla każdego wejściowego pliku bam w katalogu ~/Desktop/GSE126612/log/filter-sort-bam.
4. Po zakończeniu uruchamiania dokładnie przejrzyj pliki dziennika. Jeśli w plikach dziennika pojawi się jakikolwiek komunikat o błędzie, napraw go i spróbuj ponownie uruchomić skrypt powłoki. W przypadku problemów z rozwiązaniem kwestii, poproś o pomoc za pomocą strony z problemami (issues) projektu Easy Shells CUTnRUN w serwisie GitHub (https://github.com/JunwooLee89/Easy-Shells_CUTnRUN/issues). 
   UWAGA: Wynikowe pliki bam (wyjściowe) posortowane według nazw fragmentów będą służyć jako pliki wejściowe do utworzenia plików fragment BED oraz raw readcounts bedGraph. Pliki bam posortowane według współrzędnych będą służyć jako pliki wejściowe do generowania plików fragment BEDPE. Wszystkie pliki BED, bedGraph i BEDPE zostaną wykorzystane do wyznaczania pików (peak calling) i wizualizacji w analizie downstream. Wszystkie pliki adnotacji BED dla kanonicznych regionów chromosomowych (chr1~22, chrX, chrY oraz chrM), publicznie adnotowanych regionów czarnej listy²³ oraz regionów powtórzeń TA¹⁸ znajdują się w katalogu ~/Desktop/Easy-Shells_CUTnRUN/blacklist. W razie potrzeby użyj tego katalogu, aby dodać dodatkowe pliki czarnej listy. Użyj tego skryptu powłoki do wykonania tych samych funkcji dla innych plików bam z zmapowanymi parami odczytów, zmieniając ścieżkę i nazwy plików bam. Wpisz w terminalu samtools view --help oraz samtools sort --help, aby uzyskać więcej informacji na temat tych funkcji.

9. Konwersja zmapowanych par odczytów do plików fragment BEDPE, BED oraz plików bedGraph z surową liczbą odczytów

1. Otwórz terminal i wpisz echo $SHELL, aby sprawdzić domyślną powłokę w aktywnym terminalu. Jeśli domyślną powłoką w bieżącym terminalu jest Bash, użytkownicy powinni zobaczyć w terminalu następującą informację: /path/to/bash (lub podobny komunikat, np. /bin/bash).
2. Jeśli domyślna powłoka nie jest powłoką Bash, ustaw Bash jako powłokę domyślną, wpisując w terminalu chsh -s $(which bash). Jeśli terminal używa powłoki Bash jako domyślnej, pomiń ten krok.
3. Wpisz w terminalu ~/Desktop/Easy-Shells_CUTnRUN/scripts/Script_08_bam-to-BEDPE-BED-bedGraph.sh lub przeciągnij plik skryptu powłoki do terminala i naciśnij Enter.
   UWAGA: Ten skrypt wykona następujące czynności: (i) uruchomi funkcje macs3 filterdup oraz awk, aby przekonwertować pliki bam posortowane według współrzędnych na pliki BEDPE fragmentów, których długość jest krótsza niż 1kb, a następnie zapisze pliki BEDPE w katalogu ~/Desktop/GSE126612/BEDPE. (ii) utworzy katalog logów (~/Desktop/GSE126612/log/bam-to-BEDPE) i zapisze plik logu dla każdego pliku fragmentów zmapowanych odczytów. (iii) uruchomi funkcje bedtools bamtobed oraz awk, cut, sort, aby przekonwertować pliki bam posortowane według nazw fragmentów na pliki BED fragmentów, których długość jest krótsza niż 1 kb. (iv) utworzy jeden folder (~/Desktop/GSE126612/bam-to-bed) i zapisze w nim pliki BED fragmentów. (v) zapisze plik logu dla każdego pliku BED fragmentów zmapowanych odczytów w katalogu logów (~/Desktop/GSE126612/log/bam-to-bed). (vi) wykona funkcję bedtools genomecov w celu wygenerowania surowych plików bedGraph z liczbą odczytów, wykorzystując pliki BED fragmentów, a wyniki zapisze w jednym folderze (~/Desktop/GSE126612/bedGraph).
4. Po zakończeniu uruchamiania dokładnie sprawdź pliki logów. W przypadku wystąpienia problemów z ich rozwiązaniem poproś o pomoc za pośrednictwem strony zgłoszeń GitHub Easy Shells CUTnRUN (https://github.com/JunwooLee89/Easy-Shells_CUTnRUN/issues). 
   UWAGA: Wyjściowe surowe pliki bedGraph z liczbą odczytów zostaną użyte jako pliki wejściowe dla programu peak caller SEACR z opcją normalizacji w sekcji 12 oraz normalizacją skalowanego ułamkowego odczytu (SFRC)²² w sekcji 10. Pliki BED fragmentów posłużą jako pliki wejściowe dla normalizacji Reads Per Million mapped reads in the negative Control (SRPMC) z wykorzystaniem Spike-in^24,25 w sekcji 10. Aby przechwycić tylko krótkie fragmenty (>100 bp) dla danych CUT&RUN czynników związanych z chromatyną, należy zmienić krok filtrowania fragmentów w tym skrypcie, a następnie przejść do kroku normalizacji. Aby porównać sygnały CUT&RUN pomiędzy fragmentami krótkimi i o regularnej wielkości w obrębie tej samej próbki, normalizacja SFRC może być pomocna w celu zmniejszenia potencjalnego efektu down-samplingu spowodowanego przechwytywaniem wyłącznie krótkich fragmentów. Użyj tego skryptu powłoki, aby przeprowadzić analogiczne procesy dla innych posortowanych plików bam z sekwencjonowaniem paired-end, zmieniając ścieżkę i format nazw plików bam oraz bed.

10. Konwersja plików bedGraph z surowymi liczbami odczytów do znormalizowanych plików bedGraph i bigWig

1. Otwórz terminal i wpisz echo $SHELL, aby sprawdzić domyślną powłokę w aktywnym terminalu. Jeśli domyślną powłoką w aktualnym terminalu jest Bash, użytkownicy mogą zobaczyć w terminalu następującą informację: /path/to/bash (lub podobny komunikat, np. /bin/bash).
2. Jeśli domyślna powłoka nie jest powłoką Bash, ustaw Bash jako powłokę domyślną, wpisując w terminalu chsh -s $(which bash). Jeśli terminal używa już powłoki Bash jako domyślnej, pomiń ten krok.
3. Wpisz ~/Desktop/Easy-Shells_CUTnRUN/scripts/Script_09_normalization_SFRC.sh w terminalu lub przeciągnij plik skryptu powłoki do terminala i naciśnij enter.
   UWAGA: Ten skrypt został napisany w celu: (i) uruchomienia pętli for z funkcją awk w celu utworzenia plików bedGraph znormalizowanych metodą SFRC, wykorzystując surowe pliki bedGraph z liczbą odczytów znajdujące się w ~/Desktop/GSE126612/bedGraph. (ii) wykonania funkcji bedGraphToBigWig w celu utworzenia skompresowanego formatu (.bw) znormalizowanych plików bedGraph SFRC w ~/Desktop/GSE126612/bigWig. (iii) zapisania pliku dziennika rejestrującego czynnik normalizacji użyty do obliczeń SFRC dla każdego uruchomienia i zapisania pliku dziennika w ~/Desktop/GSE126612/log/SFRC.
4. Po zakończeniu uruchomienia sprawdź pliki dziennika. Jeśli pojawi się jakikolwiek komunikat o błędzie, popraw błąd i ponownie uruchom skrypt powłoki. W przypadku problemów z rozwiązaniem kwestii, poproś o pomoc za pomocą strony problemów (issues) projektu Easy Shells CUTnRUN na GitHubie (https://github.com/JunwooLee89/Easy-Shells_CUTnRUN/issues). 
   UWAGA: Normalizacja skalowanej ułamkowej liczby odczytów (scaled fractional readcount normalization) została zastosowana w oryginalnej publikacji²² zbioru danych CUT&RUN GSE126612. Wzór normalizacji dla bin i jest następujący:
   
   Ponieważ ta metoda normalizacji nie obejmuje normalizacji z kontrolą negatywną (na przykład próbką IgG) ani kontrolą spike-in, podejście to może nie być idealne do obserwacji różnic w sygnale w całym genomie między próbkami. Jednakże, ponieważ metoda ta jest teoretycznie podobna do innych normalizacji opartych na całkowitej liczbie odczytów (na przykład Count Per Million), powinna być wystarczająca do obserwacji lokalnych różnic w sygnale między próbkami.
5. Wpisz ~/Desktop/Easy-Shells_CUTnRUN/scripts/Script_09_normalization_SRPMC.sh w terminalu lub przeciągnij plik skryptu powłoki do terminala i naciśnij enter.
   UWAGA: Ten skrypt będzie: (i) uruchamiać pętlę for z funkcją bedtools genomecov w celu utworzenia plików bedgraph znormalizowanych metodą SRPMC w ~/Desktop/GSE126612/bedGraph, wykorzystując pliki fragment BED w ~/Desktop/GSE126612/bam-to-bed. (ii) zapisywać plik dziennika rejestrujący czynniki normalizacji użyte do normalizacji SRPMC dla każdego uruchomienia w ~/Desktop/GSE126612/log/SRPMC. (iii) wykonywać funkcję bedGraphToBigWig w celu utworzenia skompresowanego formatu (.bw) znormalizowanych plików bedGraph i zapisywania znormalizowanych plików bigWig w folderze ~/Desktop/GSE126612/bigWig.
6. Po zakończeniu uruchomienia dokładnie przejrzyj pliki dziennika. Jeśli w plikach dziennika znajduje się jakikolwiek komunikat o błędzie, popraw błąd i ponownie uruchom skrypt powłoki. W przypadku problemów z rozwiązaniem kwestii, poproś o pomoc za pomocą strony problemów (issues) projektu Easy Shells CUTnRUN na GitHubie (https://github.com/JunwooLee89/Easy-Shells_CUTnRUN/issues).
   UWAGA: Wzór normalizacji SRPMC został opracowany w celu normalizacji rzeczywistej liczby odczytów próbek zarówno z kontrolą negatywną (np. próbka IgG), jak i kontrolą spike-in, poprzez połączenie czynnika normalizacji RPM (Reads Per Million mapped reads), RPS (stosunek odczytów do odczytów spike-in) oraz względnego stosunku sygnału do kontroli^24,25. Definicja RPS jest następująca:
   
   Poprzez zastosowanie RPS zarówno dla rzeczywistej próbki, jak i próbki kontrolnej negatywnej, względny stosunek sygnału (RS) do kontroli dla rzeczywistej próbki można obliczyć następująco:
   
   A definicja czynnika normalizacji RPM (RPM:NF) jest następująca:
   
   Na tej podstawie czynnik normalizacji SRPMC (SRPMC:NF) został uzyskany poprzez połączenie RS i RPM:NF:
   
   Wzór ten można uprościć do postaci:
   
   Zatem metoda SRPMC normalizuje odczyty na podstawie (1) stosunku odczytów spike-in między kontrolą a próbką oraz (2) odczytów kontrolnych znormalizowanych RPM. Ponieważ czynnik normalizacji ten uwzględnia odczyty spike-in i sprawia, że odczyty kontrolne są porównywalne między próbkami, metoda ta byłaby odpowiednia do obserwacji różnic w całym genomie między próbkami i redukcji efektu serii (batch effect) w całkowitej liczbie odczytów rzeczywistych próbek i kontroli w różnych seriach eksperymentów. Znormalizowane pliki bedGraph staną się plikami wejściowymi do wywoływania piki za pomocą SEACR w sekcji 11. A te znormalizowane pliki bigWig będą wykorzystane w wizualizacji loci poprzez IGV oraz w tworzeniu map ciepła (heatmap) i wykresów średnich za pomocą Deeptools. Zdecydowanie sugeruje się użycie przeglądarki genomu do wizualizacji wzorca krajobrazu zbioru danych CUT&RUN przy użyciu znormalizowanych plików bigWig w reprezentatywnych regionach genomowych w celu oceny jakości danych. Próbki CUT&RUN wykazujące szumne wzorce sygnału tła przypominające kontrolę IgG powinny zostać prawdopodobnie pominięte w dalszych analizach. Użyj tych skryptów powłoki do normalizacji innych plików bed z odczytami i surowych plików bedGraph z liczbą odczytów, zmieniając ścieżki i nazwy plików dla plików bed i bedgraph wejściowych oraz wyjściowych. Edytuj te skrypty, aby zastosować inne obliczenia normalizacji, zmieniając czynniki i wzory wewnątrz skryptu.

11. Walidacja rozkładu wielkości fragmentów

1. Otwórz terminal i wpisz echo $SHELL , aby sprawdzić domyślną powłokę w aktywnym terminalu. Jeśli domyślną powłoką jest Bash, w terminalu pojawi się następujący komunikat: /path/to/bash (lub podobny, np. /bin/bash).
2. Jeśli domyślna powłoka nie jest powłoką Bash, ustaw ją jako domyślną, wpisując w terminalu chsh -s $(which bash). Jeśli terminal domyślnie używa powłoki Bash, pomiń ten krok.
3. Wpisz w terminalu ~/Desktop/Easy-Shells_CUTnRUN/scripts/Script_10_insert-size-analysis.sh lub przeciągnij plik skryptu powłoki do terminala i naciśnij Enter.
   UWAGA: Skrypt ten został napisany w celu: (i) uruchomienia funkcji picard.jar CollectInsertSizeMetrics przy użyciu plików BAM z zmapowanymi parami odczytów znajdujących się w folderze ~/Desktop/GSE126612/filtered-bam, aby określić rozkład wielkości wstawek. (ii) utworzenia jednego folderu (~/Desktop/GSE126612/insert-size-distribution) i zapisania w nim wyników analizy rozkładu wielkości wstawek. (iii) wygenerowania pliku dziennika dla każdego wejściowego pliku BAM w folderze ~/Desktop/GSE126612/log/insert-size-distribution .
4. Po zakończeniu uruchamiania dokładnie sprawdź pliki dziennika. Jeśli w plikach dziennika pojawi się jakikolwiek komunikat o błędzie, popraw go i spróbuj ponownie uruchomić skrypt powłoki. W przypadku problemów z rozwiązaniem błędu poproś o pomoc za pośrednictwem strony zgłoszeń GitHub projektu Easy Shells CUTnRUN (https://github.com/JunwooLee89/Easy-Shells_CUTnRUN/issues). 
   UWAGA: Zazwyczaj analiza wielkości wstawek (wynik) dla próbek CUT&RUN wykazuje główne piki w zakresach wielkości mononukleosomowych (100-300 bp) i dinukleosomowych (300-500 bp). Błędy techniczne/ograniczenia (takie jak nadmierna lub niewystarczająca trawienie przez MNazę podczas przygotowania próbek CUT&RUN lub niewłaściwa selekcja wielkości podczas przygotowania biblioteki) mogą powodować wzbogacenie fragmentów o wielkości równej lub większej niż trinukleosomowe (500-700 bp) oraz równej lub mniejszej niż subnukleosomowe (<100 bp). Czasami brak pików o wielkości mononukleosomowej przy jednoczesnym wzbogaceniu fragmentów długich (>500 bp) i krótkich (<100 bp) może wynikać z zakresów selekcji wielkości biblioteki wybranych na etapie pracy w laboratorium mokrym lub niskiej głębokości sekwencjonowania. Aby ocenić jakość przetworzonych próbek CUT&RUN, porównaj ze sobą głębokość sekwencjonowania („total sequenced bases” / „total reference genome size”), ogólny obraz genomowy z wykorzystaniem znormalizowanych plików bigWig z liczbą odczytów z sekcji 10 oraz wzorzec rozkładu wielkości wstawek. Linie przerywane na histogramach reprezentują „cumulative fraction” odczytów o wielkości wstawki większej lub równej wartości na osi x. Linia przerywana ta umożliwia identyfikację rozkładu wielkości wstawek w wejściowym pliku zmapowanych odczytów. Postęp wzdłuż osi x wiąże się ze zwiększającą się wielkością wstawki. Linia przerywana wskazuje proporcję zmapowanych par odczytów w wejściowym pliku BAM, które mają wielkość wstawki co najmniej tak dużą, jak wskazano w przecinającym punkcie osi x. Zatem interpretacja zaczyna się od 1 po lewej stronie, co oznacza, że wszystkie odczyty mają wielkość wstawki większą lub równą najmniejszemu rozmiarowi, i maleje w stronę 0 wraz ze wzrostem wielkości wstawki.

12. Wyznaczanie szczytów (peak calling) przy użyciu programów MACS2, MACS3 i SEACR

1. Otwórz terminal i wpisz echo $SHELL, aby sprawdzić domyślną powłokę w aktywnym terminalu. Jeśli domyślną powłoką w aktualnym terminalu jest Bash, użytkownicy mogą zobaczyć w terminalu następującą ścieżkę: /path/to/bash (lub podobny komunikat, np. /bin/bash).
2. Jeśli domyślna powłoka nie jest powłoką Bash, ustaw Bash jako powłokę domyślną, wpisując w terminalu chsh -s $(which bash). Jeśli terminal używa powłoki Bash jako domyślnej, pomiń ten krok.
3. Wpisz w terminalu ~/Desktop/Easy-Shells_CUTnRUN/scripts/Script_11_peak-calling_MACS.sh lub przeciągnij plik skryptu powłoki do terminala i naciśnij enter.
   UWAGA: Ten skrypt został napisany w celu: (i) uruchomienia funkcji macs2 callpeak oraz macs3 callpeak z kontrolą IgG i bez niej, wykorzystując pliki BEDPE fragmentów do wyznaczania pików (peak calling) oraz zapisania wyników wyznaczania pików w katalogach wyjściowych (~/Desktop/GSE126612/MACS2 oraz ~/Desktop/GSE126612/MACS3). (ii) zapisania logów tych procesów wyznaczania pików w formie plików tekstowych w katalogu logów (~/Desktop/GSE126612/log/MACS2 oraz ~/Desktop/GSE126612/log/MACS3)
4. Wpisz w terminalu ~/Desktop/Easy-Shells_CUTnRUN/scripts/Script_11_peak-calling_SEACR.sh lub przeciągnij plik skryptu powłoki do terminala i naciśnij enter.
   UWAGA: Ten skrypt został napisany w celu: (i) uruchomienia skryptu SEACR_1.3.sh z kontrolą IgG i bez niej, z opcjami rygorystyczną (stringent) i łagodną (relaxed), wykorzystując pliki bedGraph z surową liczbą odczytów (raw readcounts) oraz znormalizowane pliki bedGraph do wyznaczania pików. (ii) utworzenia katalogu wyjściowego (~/Desktop/GSE126612/SEACR-peaks) i zapisania wyników wyznaczania pików przez SEACR. (iii) zapisania logów tych procesów wyznaczania pików w formie plików tekstowych w katalogu logów (~/Desktop/GSE126612/log/SEACR).
5. Po zakończeniu uruchamiania skryptów powłoki dokładnie sprawdź pliki logów. Jeśli w plikach logów pojawi się jakikolwiek komunikat o błędzie, najpierw napraw błąd. Niektóre programy mogą nie wyznaczać pików dla próbki kontrolnej IgG przy jednoczesnym zastosowaniu opcji kontroli IgG, dlatego należy zignorować komunikaty o błędach dotyczące próbki kontrolnej IgG z włączoną opcją kontroli IgG. W przypadku problemów z rozwiązaniem kwestii, poproś o pomoc za pośrednictwem strony z problemami (issues) projektu Easy Shells CUTnRUN na GitHubie (https://github.com/JunwooLee89/Easy-Shells_CUTnRUN/issues). 
   UWAGA: Te dwa skrypty powłoki wykonują wyznaczanie pików dla próbek CUT&RUN przy użyciu trzech programów (MACS2, MACS3 i SEACR) z różnymi opcjami: z opcją kontroli IgG lub bez niej, przy użyciu plików bedGraph z surową liczbą odczytów z opcją normalizacji programu do wyznaczania pików lub znormalizowanych plików bedGraph bez opcji normalizacji, a także z rygorystycznymi i łagodnymi opcjami wyznaczania pików w SEACR. Ponieważ pliki wyjściowe z wyznaczania pików nie nadają się do bezpośredniego użycia w analizach downstream, Easy Shells CUTnRUN zawiera jeden skrypt do przetwarzania tych plików wyjściowych w celu utworzenia nowych plików pików, które zawierają chromosom, początek, koniec i nazwę pików. Dzięki zastosowaniu intensywnych podejść do wyznaczania pików, Easy Shells CUTnRUN daje możliwość wyboru programu najlepiej dopasowanego do projektu CUT&RUN użytkownika poprzez porównanie pików wyznaczonych przez trzy różne programy. Dodatkowo, ten rurociąg analityczny CUT&RUN umożliwia wybór opcji wyznaczania pików najlepiej dopasowanych do projektu CUT&RUN użytkownika. Porównania te zostaną przeprowadzone za pomocą diagramów Venna oraz wizualizacji w postaci map ciepła (heatmap) i wykresów średnich.

13. Tworzenie plików BED dla zidentyfikowanych pików

1. Otwórz terminal i wpisz echo $SHELL, aby sprawdzić domyślną powłokę w aktywnym terminalu. Jeśli domyślną powłoką w bieżącym terminalu jest Bash, użytkownicy mogą zobaczyć w terminalu następującą informację: /path/to/bash (lub podobny komunikat, taki jak /bin/bash).
2. Jeśli domyślna powłoka nie jest powłoką Bash, ustaw Bash jako powłokę domyślną, wpisując w terminalu chsh -s $(which bash). Jeśli terminal używa powłoki Bash jako domyślnej, pomiń ten krok.
3. Wpisz ~/Desktop/Easy-Shells_CUTnRUN/scripts/Script_12_make-peak-bed.sh w terminalu lub przeciągnij plik skryptu powłoki do terminala i naciśnij enter.
   UWAGA: Ten skrypt służy do: (i) uruchomienia funkcji awk z wykorzystaniem plików bed w folderze ~/Desktop/GSE126612/SEACR w celu utworzenia dwóch typów plików bed dla pików SEACR w folderze ~/Desktop/GSE126612/peak-bed_SEACR. Pełne pliki bed dla pików zawierają początek i koniec każdego piku, a skoncentrowane pliki bed dla pików zawierają początek i koniec bin z najwyższym sygnałem w obrębie każdego piku. (ii) uruchomienia funkcji awk z wykorzystaniem plików ***_peaks.xls w folderach ~/Desktop/GSE126612/MACS2 oraz ~/Desktop/GSE126612/MACS3 w celu utworzenia pełnych plików bed dla pików, które zawierają początek i koniec każdego piku wykrytego przez MACS2 i MACS3, w folderach ~/Desktop/GSE126612/peak-bed_MACS2 oraz ~/Desktop/GSE126612/peak-bed_MACS3 . (iii) uruchomienia funkcji awk z wykorzystaniem plików ***_summits.bed w folderach ~/Desktop/GSE126612/MACS2 oraz ~/Desktop/GSE126612/MACS3 w celu utworzenia skoncentrowanych plików bed dla pików, które zawierają początek i koniec najbardziej znaczącego binu w obrębie każdego piku. (iv) pliki dziennika są zapisywane w formacie tekstowym w folderze ~/Desktop/GSE126612/log/peak-bed .
4. Wpisz ~/Desktop/Easy-Shells_CUTnRUN/scripts/Script_13_filter-peaks.sh w terminalu lub przeciągnij plik skryptu powłoki do terminala i naciśnij enter.
   UWAGA: Ten skrypt służy do: (i) uruchomienia funkcji bedtools intersect z wykorzystaniem plików bed dla pików, które zostały wywołane bez opcji kontroli IgG, aby usunąć piki nakładające się na piki kontroli IgG. (ii) przefiltrowane pliki bed dla pików są zapisywane w folderach ~/Desktop/GSE126612/peak-bed-filtered_MACS2, ~/Desktop/GSE126612/peak-bed-filtered_MACS3 oraz ~/Desktop/GSE126612/peak-bed-filtered_SEACR . (iii) plik dziennika log_filter-peaks.txt jest tworzony w folderze ~/Desktop/GSE126612/log/filter-peaks.
5. Wpisz ~/Desktop/Easy-Shells_CUTnRUN/scripts/Script_14_cat-merge-peak-bed_MACS.sh w terminalu lub przeciągnij plik skryptu powłoki do terminala i naciśnij enter.
   UWAGA: Ten skrypt służy do: (i) uruchomienia funkcji cat oraz sort w celu połączenia pełnych plików bed dla pików z powtórzeń MACS2 i MACS3 w jeden plik bed dla pików oraz posortowania połączonego pliku w folderze ~/Desktop/GSE126612/bed-for-comparison. (ii) uruchomienia funkcji bedtools merge z wykorzystaniem połączonych pełnych plików bed dla pików w celu scalenia pików, które na siebie nachodzą. (iii) plik dziennika log_cat-merged-peak-bed_MACS.txt jest zapisywany w folderze dzienników ~/Desktop/GSE126612/log/cat-merged-peak-bed.
6. Wpisz ~/Desktop/Easy-Shells_CUTnRUN/scripts/Script_14_cat-merge-peak-bed_SEACR.sh w terminalu lub przeciągnij plik skryptu powłoki do terminala i naciśnij enter.
   UWAGA: Ten skrypt służy do: (i) uruchomienia funkcji cat oraz sort w celu połączenia pełnych plików bed dla pików z powtórzeń SEACR w jeden plik bed dla pików oraz posortowania połączonego pliku w folderze ~/Desktop/GSE126612/bed-for-comparison . (ii) uruchomienia funkcji bedtools merge z wykorzystaniem połączonych pełnych plików bed dla pików w celu scalenia pików, które na siebie nachodzą. (iii) plik dziennika log_cat-merged-peak-bed_SEACR.txt jest zapisywany w folderze dzienników ~/Desktop/GSE126612/log/cat-merged-peak-bed.
7. Po zakończeniu uruchamiania skryptów powłoki dokładnie przejrzyj pliki dziennika. Jeśli w plikach dziennika pojawi się jakikolwiek komunikat o błędzie, napraw błąd i uruchom skrypt(y) ponownie. W przypadku problemów z rozwiązaniem kwestii, poproś o pomoc za pośrednictwem strony zgłaszania problemów projektu Easy Shells CUTnRUN na GitHubie (https://github.com/JunwooLee89/Easy-Shells_CUTnRUN/issues).
   UWAGA: Pełne pliki bed dla obszarów pików zostaną użyte jako pliki wejściowe do analizy diagramu Venna w celu porównania podobieństwa między opcjami wykrywania pików, metodami wykrywania pików, powtórzeniami oraz obserwacjami krajobrazu genomicznego w pobliżu obszarów pików. Scalone pełne pliki bed dla obszarów pików zostaną wykorzystane do analizy głównych składowych (PC) i analizy korelacji współczynnika Pearsona przy użyciu deeptools. Skoncentrowane pliki bed dla pików zostaną wykorzystane do analizy map ciepła (heatmap) i wykresów średnich przy użyciu Deeptools.

14. Walidacja podobieństwa między powtórzeniami z wykorzystaniem korelacji Pearsona i analizy głównych składowych (PC).

1. Otwórz terminal i wpisz echo $SHELL, aby sprawdzić domyślną powłokę w aktywnym terminalu. Jeśli domyślną powłoką w bieżącym terminalu jest Bash, użytkownicy mogą zobaczyć w terminalu komunikat: /path/to/bash (lub podobny, np. /bin/bash).
2. Jeśli domyślna powłoka nie jest powłoką Bash, ustaw Bash jako powłokę domyślną, wpisując w terminalu chsh -s $(which bash). Jeśli terminal korzysta już z powłoki Bash jako domyślnej, pomiń ten krok.
3. Wpisz w terminalu ~/Desktop/Easy-Shells_CUTnRUN/scripts/Script_15_correlation_plotCorrelation.sh lub przeciągnij plik skryptu powłoki do terminala i zatwierdź klawiszem enter.
   UWAGA: Skrypt ten został napisany w celu: (i) uruchomienia funkcji multiBamSummary BED-file z wykorzystaniem plików bam z powtórzeń, które zostały posortowane według współrzędnych, oraz połączonych plików bed wszystkich szczytów dla CTCF, H3K27Ac i RNAPII-S5P, aby wygenerować pliki macierzy do analizy korelacji Pearsona w folderze Desktop/GSE126612/deeptools_multiBamSummary. (ii) uruchomienia funkcji plotCorrelation z wykorzystaniem plików macierzy w celu obliczenia współczynnika korelacji Pearsona i przeprowadzenia klastrowania na mapie ciepła (heatmap), a następnie zapisania wyniku w folderze ~/Desktop/GSE126612/deeptools_plotCorrelation. (iii) zapisania pliku dziennika log_plotCorrelation.txt w folderze ~/Desktop/GSE126612/log/correlation.
4. Wpisz w terminalu ~/Desktop/Easy-Shells_CUTnRUN/scripts/Script_15_correlation_plotPCA.sh lub przeciągnij plik skryptu powłoki do terminala i zatwierdź klawiszem enter.
   UWAGA: Skrypt ten został napisany w celu: (i) uruchomienia funkcji multiBamSummary BED-file z wykorzystaniem plików bam, które zostały posortowane według współrzędnych, oraz połączonych plików bed wszystkich szczytów, obejmujących wszystkie szczyty CTCF, H3K27ac i RNAPII-S5P, aby wygenerować pliki macierzy do analizy głównych składowych (PCA) w folderze Desktop/GSE126612/deeptools_multiBamSummary. (ii) uruchomienia funkcji plotPCA z wykorzystaniem plików macierzy w celu przeprowadzenia PCA i zapisania wyniku w folderze ~/Desktop/GSE126612/deeptools_plotPCA. (iii) zapisania pliku dziennika log_plotPCA.txt w folderze ~/Desktop/GSE126612/log/correlation.
5. Po zakończeniu działania skryptów powłoki sprawdź pliki dziennika. W przypadku wystąpienia jakiegokolwiek komunikatu o błędzie, popraw błąd i uruchom skrypty powłoki ponownie. Jeśli pojawią się problemy z rozwiązaniem problemu, poproś o pomoc za pośrednictwem strony z problemami (issues) w serwisie github dla Easy Shells CUTnRUN (https://github.com/JunwooLee89/Easy-Shells_CUTnRUN/issues). 
   UWAGA: Co do zasady, prawidłowo przygotowane i przetworzone powtórzenia wykazują wyższe wartości współczynnika korelacji Pearsona w obrębie tej samej grupy klastrowania oraz bliskie rozmieszczenie w analizie głównych składowych. Każde powtórzenie, które wykazuje niższy współczynnik korelacji Pearsona i dużą odległość od innych powtórzeń na wykresie głównych składowych, może reprezentować potencjalną wartość odstającą wśród powtórzeń. Ten skrypt powłoki jest odpowiedni dla wszelkich danych z odczytami zmapowanymi w formacie bam. Zmień ścieżkę i nazwę plików bigwig, aby spełnić specyficzne wymagania projektu.

15. Walidacja podobieństwa między powtórzeniami, metodami i opcjami wyznaczania pików za pomocą diagramu Venna

1. Otwórz terminal i wpisz echo $SHELL, aby sprawdzić domyślną powłokę w aktywnym terminalu. Jeśli domyślną powłoką w bieżącym terminalu jest Bash, w terminalu pojawi się ścieżka taka jak /path/to/bash (na przykład /bin/bash).
2. Jeśli domyślna powłoka to nie Bash, ustaw powłokę Bash jako domyślną, wpisując w terminalu chsh -s $(which bash). Jeśli terminal używa powłoki Bash jako domyślnej, pomiń ten krok
3. Wpisz w terminalu ~/Desktop/Easy-Shells_CUTnRUN/scripts/Script_16_venn-diagram_methods.sh lub przeciągnij plik skryptu powłoki do terminala i naciśnij enter.
   UWAGA: Ten skrypt został napisany, aby: (i) uruchomić funkcję intervene venn przy użyciu plików bed z całymi regionami peaków, aby znaleźć nakładanie się peaków wywołanych przez różne opcje (z opcją kontroli IgG lub bez niej, z normalizacją lub bez niej oraz przy rygorystycznych lub łagodnych opcjach wywołania peaków dla SEACR). (ii) utworzyć jeden folder (~/Desktop/GSE126612/intervene_methods) i zapisać w nim wyniki analizy diagramu Venna. (iii) utworzyć plik logu log_intervene_methods.txt w folderze ~/Desktop/GSE126612/log/intervene.
4. Wpisz w terminalu ~/Desktop/Easy-Shells_CUTnRUN/scripts/Script_16_venn-diagram_replicates.sh lub przeciągnij plik skryptu powłoki do terminala i naciśnij enter.
   UWAGA: Ten skrypt został napisany, aby: (i) uruchomić funkcję intervene venn przy użyciu plików bed z całymi regionami peaków, aby znaleźć nakładanie się peaków pomiędzy replikatami. (ii) utworzyć jeden folder (~/Desktop/GSE126612/intervene_replicates) i zapisać w nim wyniki analizy diagramu Venna. (iii) utworzyć plik logu log_intervene_replicates.txt w folderze ~/Desktop/GSE126612/log/intervene.
5. Po zakończeniu uruchamiania skryptów powłoki przejrzyj pliki logów. W przypadku wystąpienia jakiegokolwiek komunikatu o błędzie, napraw błąd i ponownie uruchom skrypty powłoki. W przypadku jakichkolwiek problemów z użyciem potoku analizy Easy Shells CUTnRUN, poproś o pomoc na stronie problemów GitHub Easy Shells CUTnRUN (https://github.com/JunwooLee89/Easy-Shells_CUTnRUN/issues).
   UWAGA: Wyniki analizy diagramu Venna pozwalają dobrać najodpowiedniejsze opcje wywoływania peaków, metody oraz replikaty o wysokiej powtarzalności do dalszych analiz. Zaleca się wybór opcji i metod wywoływania peaków, które wykazują najwyższą liczbę wykrytych peaków przy dobrym nakładaniu się z innymi metodami i opcjami wywoływania peaków.

16. Analiza map ciepła i wykresów średnich w celu wizualizacji zidentyfikowanych pików.

1. Otwórz terminal i wpisz echo $SHELL , aby sprawdzić domyślną powłokę w aktywnym terminalu. Jeśli domyślną powłoką w bieżącym terminalu jest Bash, w terminalu pojawi się ścieżka taka jak /path/to/bash (na przykład /bin/bash).
2. Jeśli domyślna powłoka nie jest powłoką Bash, ustaw Bash jako powłokę domyślną, wpisując w terminalu chsh -s $(which bash). Jeśli terminal korzysta z powłoki Bash jako domyślnej, można pominąć ten krok.
3. Wpisz w terminalu ~/Desktop/Easy-Shells_CUTnRUN/scripts/Script_27_plotHeatmap_focused.sh lub przeciągnij plik skryptu powłoki do terminala i naciśnij enter.
   UWAGA: Skrypt ten służy do: (i) uruchomienia funkcji computeMatrix reference-point z wykorzystaniem znormalizowanych plików bigWig oraz plików bed dla wybranych szczytów (focused peaks), aby utworzyć znormalizowane macierze liczby odczytów w centrum wybranych szczytów w folderze ~/Desktop/GSE126612/deeptools_computeMatrix. (ii) uruchomienia funkcji plotHeatmap z wykorzystaniem znormalizowanej macierzy liczby odczytów w celu wygenerowania map ciepła i wykresów średnich, które wizualizują wzorzec rozkładu znormalizowanej liczby odczytów w lokalizacjach wybranych szczytów. (iii) utworzenia jednego folderu (~/Desktop/GSE126612/deeptools_plotHeatmap) i zapisania w nim plików wyjściowych funkcji plotHeatmap. (iv) zapisania jednego pliku dziennika log_plotHeatmap_focused.txt w folderze ~/Desktop/GSE126612/log/plotHeatmap.
4. Wpisz w terminalu ~/Desktop/Easy-Shells_CUTnRUN/scripts/Script_27_plotHeatmap_whole.sh lub przeciągnij plik skryptu powłoki do terminala i naciśnij enter.
   Skrypt ten służy do: (i) uruchomienia funkcji computeMatrix reference-point z wykorzystaniem znormalizowanych plików bigWig oraz plików bed dla wszystkich szczytów (whole peaks), aby utworzyć znormalizowane macierze liczby odczytów w centrum wszystkich szczytów w folderze ~/Desktop/GSE126612/deeptools_computeMatrix. (ii) uruchomienia funkcji plotHeatmap z wykorzystaniem znormalizowanej macierzy liczby odczytów w celu wygenerowania map ciepła i wykresów średnich, które wizualizują wzorzec rozkładu znormalizowanej liczby odczytów w lokalizacjach wszystkich szczytów. (iii) utworzenia jednego folderu (~/Desktop/GSE126612/deeptools_plotHeatmap) i zapisania w nim plików wyjściowych funkcji plotHeatmap. (iv) zapisania jednego pliku dziennika log_plotHeatmap_whole.txt w folderze ~/Desktop/GSE126612/log/plotHeatmap.
5. Po zakończeniu działania skryptów powłoki przejrzyj pliki dziennika. W przypadku wystąpienia jakiegokolwiek komunikatu o błędzie, napraw go i ponownie uruchom skrypty. W razie problemów z korzystaniem z potoku analizy Easy Shells CUTnRUN, poproś o pomoc na stronie zgłaszania błędów projektu Easy Shells CUTnRUN w serwisie GitHub (https://github.com/JunwooLee89/Easy-Shells_CUTnRUN/issues).
   UWAGA: W idealnych warunkach lokalizacje wierzchołków szczytów MACS2/3 oraz lokalizacje wybranych szczytów SEACR wykazują ostry i skoncentrowany rozkład sygnału w centrum wykresów. Jeśli jednak algorytm wyznaczania szczytów nie działa prawidłowo dla danych CUT&RUN, na wykresach może pojawić się mniej skoncentrowany, „zaszumiony” rozkład sygnału. Zatem liczba wyznaczonych szczytów oraz wzorce rozkładu sygnału szczytów na wykresach wyjściowych będą pomocne przy ocenie poprawności szczytów do dalszej analizy CUT&RUN, obejmującej m.in. adnotację szczytów.

Wyniki

Kontrola jakości i przycinanie adapterów pozwala na zachowanie odczytów o wysokiej jakości sekwencjonowania
Techniki sekwencjonowania wysokoprzepustowego są podatne na powstawanie błędów sekwencjonowania, takich jak „mutacje” sekwencji w odczytach. Ponadto w zbiorach danych sekwencyjnych mogą występować wzbogacone dimery adapterów sekwencjonujących z powodu niewystarczającego usunięcia adapterów podczas przygotowywania biblioteki. Nadmierne błędy sekwencjonowania, takie jak mutacje odczytów, generowanie odczytów krótszych niż wymagane do prawidłowego mapowania oraz wzbogacenie dimerów adapterów, mogą wydłużyć czas mapowania odczytów i prowadzić do powstania fałszywie dodatnich zmapowanych odczytów, co zniekształca wyniki późniejszych analiz bioinformatycznych. Dlatego niezbędne jest filtrowanie jakościowe i przycinanie adapterów, aby zachować odczyty wysokiej jakości do dalszych analiz i interpretacji.

Aby zachować wysoką jakość odczytów do analizy, ten potok analizy CUT&RUN (Rysunek 2) wykorzystuje narzędzia FastQC²⁶ oraz Trim Galore²⁷. Skrypt powłoki „Script_03_fastQC.sh” uruchamia FastQC dla wszystkich plików fastq w katalogu roboczym. Wyniki (Rysunek 3) tego etapu, uzyskane przy użyciu publicznie dostępnego zbioru danych CTCF CUT&RUN z GSE126612 (SRR8581589), wskazują na obecność odczytów z zasadami o niskiej jakości (Rysunek 3A,C) oraz pewien stopień niezgodności w rozkładzie zawartości GC na sekwencję pomiędzy szacunkiem teoretycznym a rzeczywistymi odczytami (Rysunek 3E).

Wykonanie skryptu "Script_04_trimming.sh" w celu uruchomienia programu Trim Galore skutecznie usuwa odczyty z zasadami o niskiej jakości (poniżej 20 na Rysunku 3A) oraz odczyty o niskiej średniej jakości sekwencji widocznej przed przycinaniem (Rysunek 3B-D). Ponadto "Script_04_trimming.sh" skutecznie usuwa wzbogacenie średniej zawartości GC na poziomie 55~60%, przedstawione na wykresie rozkładu GC w sekwencji przed przycinaniem (Rysunek 3E,F). Wyniki te dowodzą, że ten potok analizy CUT&RUN filtruje odczyty wysokiej jakości, aby ułatwić szybkie i dokładne mapowanie odczytów do genomu referencyjnego.

Rozkład wielkości wstawek może dostarczyć szacunkowych danych dla wyników wyznaczania pików (peak calling)
Ze względu na zastosowanie MNazy w metodzie CUT&RUN (Rysunek 1), oczekuje się, że zmapowane odczyty CUT&RUN będą wykazywać piki wielkości fragmentów DNA o charakterze mononukleosomalnym (~200 bp) i dinukleosomalnym (~350 bp) na wykresach rozkładu wielkości wstawek (Rysunek 4). Problemy z detekcją niektórych celów mogą skutkować krótkimi wstawkami (< 100 bp) (Rysunek 4C). Wysoki poziom krótkich odczytów zmniejsza liczbę odczytów, które mogą zostać wykorzystane do wysokopoziomowego wyznaczania pików, co redukuje liczbę pików i wpływa na dalszą analizę. W tym potoku analizy CUT&RUN skrypt „Script_10_insert-size-analysis.sh” uruchamia funkcję „picard.jar CollectInsertSizeMetrics”, aby przeprowadzić analizę rozkładu wielkości wstawek i wyeksportować histogramy jako wynik wizualizacji (Rysunek 2). Na wykresach wyjściowych (Rysunek 4A-C) oś X przedstawia zakres wielkości wstawek, lewa strona osi Y oraz wypełniony histogram reprezentują liczbę wstawek o wartości wskazanej na osi X, a prawa strona osi Y oraz linia przerywana wskazują skumulowaną frakcję wstawek o wielkości równej lub większej od wartości na osi X. Zatem zarówno miejsce na osi X z najbardziej gwałtowną zmianą nachylenia linii przerywanej, jak i najwyższy poziom histogramu identyfikują główną wielkość wstawki w próbce. Wśród odczytów zmapowanych na genom referencyjny (ludzki, hg19), fragmenty próbki H3K27Ac (aktywna marka histonowa) wykazują oczekiwany rozkład wielkości wstawek CUT&RUN z najwyższym pikiem o wielkości mononukleosomalnej i wykrywalnym pikiem o wielkości dinukleosomalnej (Rysunek 4B). Fragmenty próbki CTCF wykazały dodatkowe grupy w regionach długości fragmentów 100~200 bp (Rysunek 4A). Podsumowując, potok analizy CUT&RUN dostarcza łatwe w użyciu skrypty powłoki do przeprowadzania analizy rozkładu wielkości wstawek po zmapowaniu odczytów na genomach referencyjnych. Analizy te stają się istotne przy szacowaniu wydajności wyznaczania pików przed dalszą analizą.

Potok analizy Easy Shells CUTnRUN zapewnia opcje filtrowania i normalizacji w celu uzyskania wiarygodnych liczb odczytów
Jednym z krytycznych punktów analizy CUT&RUN jest uzyskanie odpowiednich zmapowanych par odczytów poprzez odfiltrowanie problematycznych par odczytów z początkowych wyników mapowania oraz normalizację przefiltrowanych liczb zmapowanych odczytów za pomocą specyficznej metody obliczeniowej, która odpowiada celom lub potrzebom analizy użytkownika. Potok analizy CUT&RUN omówiony w niniejszym badaniu zawiera skrypt „Script_07_filter-sort-bam.sh”, służący do usuwania par odczytów zmapowanych na chromosomy niekanoniczne, publicznie adnotowane regiony czarnej listy²³ oraz regiony powtórzeń TA^18,22 z par odczytów, które zostały zmapowane programem bowtie2 przy użyciu skryptu „Script_06_bowtie2-mapping.sh”. Filtrowanie to jest niezbędne do usunięcia par odczytów, które mogą generować fałszywie dodatnie wyniki, sygnały odstające (spike signals) oraz błędnie zidentyfikowane piki w dalszej analizie (Rysunek 5; regiony w żółtych ramkach).

Oprócz filtracji, zastosowanie odpowiedniej metody normalizacji jest istotnym czynnikiem umożliwiającym dokładną wizualizację różnic w sygnale pomiędzy próbkami. Dlatego potok analizy CUT&RUN zawiera skrypty „Script_09_normalization_SFRC.sh” oraz „Script_09_normalization_SRPMC.sh”, które zapewniają dwie publicznie zweryfikowane metody normalizacji – skalowaną frakcyjną liczbę odczytów (SFRC)²² oraz normalizację Spike-in w odniesieniu do liczby odczytów na milion zmapowanych odczytów w kontroli negatywnej (SRPMC)^24,25 (Rysunek 5A-D). Ponieważ SFRC nie uwzględnia w wzorze próbki kontrolnej (np. IgG) ani próbki spike-in, normalizacja SFRC może być stosowana do próbek, które nie zawierają żadnej próbki kontrolnej lub w których oczekuje się różnic w sygnale jedynie w regionach lokalnych, bez różnic w skali całego genomu. Próbki znormalizowane metodą SFRC przetworzone przez potok analizy CUT&RUN (Rysunek 5A-D; ścieżki czerwone) wykazują takie same wzorce rozkładu sygnału jak publicznie dostępne zmapowane odczyty z GEO (Rysunek 5A-D; ścieżki czarne), co sugeruje, że ten potok analizy może odtworzyć wyniki publikacyjne.

Metoda SRPMC jest przydatna do normalizacji próbek, które obejmują zarówno próbki kontrolne, jak i próbki z dodatkiem wewnętrznym (spike-in) i w których spodziewana jest globalna różnica sygnału między próbkami (Rysunek 5A-D; ścieżki zielone). Ponieważ jedna próbka H3K27Ac (SRR8581599) wykazuje znacznie wyższy stosunek „(rzeczywiste odczyty CUT&RUN)/(odczyty spike-in)” (RPS próbki; 997) niż pozostałe powtórzenia (237, 175 i 161), względne sygnały H3K27Ac wydają się różne między powtórzeniami w próbkach normalizowanych metodami SFRC i SRPMC (Rysunek 5A-D; porównanie H3K27Ac we wszystkich ścieżkach). Próbki RNAPII-S5P wykazują stosunkowo niższy RPS próbki (1,7, 0,8, 2,1) niż kontrola IgG (259), w związku z czym po normalizacji SRPMC próbki RNAPII-S5P wykazują niższy sygnał niż kontrola IgG (Rysunek 5A-D; porównanie RNAPII-S5P we wszystkich ścieżkach). Dlatego omówiony tutaj schemat analizy CUT&RUN zaleca stosowanie metody SRPMC wyłącznie dla próbek, które posiadają wystarczającą liczbę odczytów w próbkach eksperymentalnych w stosunku do odczytów kontroli IgG oraz kontroli spike-in.

Porównanie za pomocą diagramu Venna może pomóc w wyborze lepszej metody i opcji wyznaczania szczytów (peak calling)
Wiele programów do wyznaczania szczytów umożliwia identyfikację istotnie wzbogaconego zajęcia białek w całym genomie. Do programów stosowanych w analizie CUT&RUN należą dotychczas przede wszystkim programy z rodziny MACS² oraz SEACR⁴. Jednakże identyfikacja najodpowiedniejszej metody i opcji wyznaczania szczytów dla danego projektu CUT&RUN może być wyzwaniem, zwłaszcza dla osób początkujących w bioinformatyce. Dlatego potok analizy CUT&RUN obejmuje etapy analizy za pomocą diagramów Venna, aby umożliwić użytkownikom porównanie podobieństw i różnic w wynikach wyznaczania szczytów pomiędzy różnymi opcjami (Script_17_intervene-options) oraz programami do wyznaczania szczytów (Script_19_intervene_methods.sh) (Rysunek 6A-H).

Z porównania zmergeowanych pików CTCF, H3K27ac oraz RNAPII-S5P, które zostały wyznaczone z opcją kontroli IgG oraz bez niej podczas etapu wyznaczania pików, wynika, że MACS2 i MACS3 wyznaczyły więcej pików z zastosowaniem opcji kontroli IgG (Rycina 6A), natomiast SEACR wyznaczył więcej pików bez opcji kontroli IgG zarówno w ustawieniach rygorystycznych, jak i rozluźnionych (Rycina 6B-D). W związku z tym potok analizy CUT&RUN sugeruje (1) zastosowanie opcji kontroli IgG dla MACS2 i MACS3, (2) wyznaczanie pików osobno dla próbek eksperymentalnych CUT&RUN oraz próbek kontrolnych IgG, a następnie odfiltrowanie pików IgG w późniejszym etapie dla programu SEACR. Porównując MACS2 i MACS3, program MACS3 wyznaczył nieco więcej pików (Rycina 6A).

Ponadto porównanie szczytów wyznaczonych przez MACS2 i MACS3 z opcją kontroli IgG oraz SEACR bez opcji kontroli IgG pokazuje, że szczyty SEACR wyznaczone z opcją rygorystyczną (stringent) pokrywają się z szczytami MACS 2 i MACS3 w większym stopniu niż szczyty SEACR wyznaczone z opcją liberalną (relaxed) (Rysunek 6E,F). Zatem wyniki potoku analizy CUT&RUN sugerują, że opcja rygorystyczna maksymalizuje spójność SEACR z wyznaczaniem szczytów w MACS. Wreszcie, diagram Venna służący do porównania pokrycia szczytów wyznaczonych przez SEACR z normalizacją dla plików bedGraph CUT&RUN z surową liczbą odczytów oraz bez normalizacji dla plików bedGraph CUT&RUN z znormalizowaną liczbą odczytów nie wykazuje różnic między metodami SFRC i SRPMC dla SEACR z opcją rygorystyczną. Szczyty SFRC charakteryzują się znacznie większą liczbą szczytów i lepszym pokryciem z szczytami z opcją normalizacji ('norm' na Rysunku 6) niż szczyty SRPMC dla SEACR z opcjami liberalnymi (Rysunek 6G,H).

Porównania statystyczne między powtórzeniami a próbkami
Wyciągnięcie dokładnych wniosków z wielu powtórzeń wymaga oceny podobieństwa tych powtórzeń. Zastosowany tutaj potok analizy CUT&RUN wykorzystuje obliczenia współczynnika korelacji statystycznej oparte na Deeptools2¹⁵, grupowanie na mapach ciepła (heatmap clustering) oraz analizę głównych składowych (PCA), aby ułatwić identyfikację próbek i powtórzeń odpowiednich do prawidłowej analizy w dalszych etapach. Grupowanie na mapach ciepła oparte na współczynniku korelacji Pearsona wykazało statystycznie istotną korelację między powtórzeniami dla CTCF, H3K27Ac i RNAPII-S5P w obrębie zidentyfikowanych regionów szczytowych (Rysunek 7A-C). Jednak analiza PCA wykazała, że jedna próbka CTCF (SRR8581590) oraz H3K27Ac (SRR8581608) znajduje się stosunkowo daleko od pozostałych powtórzeń (Rysunek 7D) we wszystkich zidentyfikowanych regionach szczytowych dla CTCF, H3K27Ac i RNAPII-S5P.

Zgodnie z diagramem Venna służącym do porównania szczytów między powtórzeniami, szczyty CTCF (SRR8581590) wykazywały najmniejszą część wspólną z innymi powtórzeniami we wszystkich trzech wynikach programów do wyznaczania szczytów (Rysunek 7E-G), natomiast szczyty H3K27Ac (SRR8581608) wykazywały najmniejszą część wspólną z innymi powtórzeniami w wynikach wyznaczania szczytów metodą SEACR (Rysunek 7F). Szczyty H3K27Ac (SRR8581608) nie wykazywały minimalnej części wspólnej z innymi powtórzeniami w wynikach wyznaczania szczytów MACS2 i MACS3 (Rysunek 7F), co może sugerować, że odległość między powtórzeniami w analizie PCA nie jest wystarczająca do zdefiniowania próbki odstającej. W związku z tym potok analizy CUT&RUN proponuje zdefiniowanie powtórzenia odstającego jako „próbki, która wykazuje niski współczynnik korelacji Pearsona w grupie klastrowania na mapie ciepła, dużą odległość od innych powtórzeń na wykresie PCA oraz najmniejszą część wspólną szczytów między powtórzeniami”.

Wyznaczanie szczytów (peak calling) ułatwia wizualizację i interpretację danych CUT&RUN
Procedura analizy CUT&RUN szczegółowo opisana w niniejszym badaniu wykorzystuje dwa rodzaje publicznie dostępnych programów do wyznaczania szczytów: rodzinę MACS oraz SEACR. Aby zoptymalizować wizualizację wyznaczonych szczytów, procedura ta wybiera bin o najwyższym sygnale jako centrum szczytu dla analiz map ciepła (heatmap) i metaplotów. Wszystkie szczyty CTCF, H3K27Ac i RNAPII-S5P wyznaczone przez programy MACS3 i SEACR wykazały bardziej ostry wzorzec rozkładu szczytów w centrum binów o najwyższym sygnale (Rycina 8A-F, wykresy „focused”) niż w centrum całych regionów szczytowych (Rycina 8A-F, wykresy „whole”). Próbki CUT&RUN przetworzone za pomocą procedury analizy Easy Shells CUTnRUN z normalizacją SFRC (Rycina 8 A-F, wykresy „SFRC”) wykazują podobne wzorce rozkładu sygnału jak próbki z normalizacją SFRC, których surowe zmapowane pary odczytów są publicznie dostępne w bazie GEO (Rycina 8A-F, wykresy „public”) w obrębie szczytów wyznaczonych przez procedurę analizy. Zatem procedura analizy CUT&RUN potrafi skutecznie odtworzyć wyniki publikowane w literaturze.

Rysunek 1: Schemat procedury eksperymentalnej CUT&RUN. CUT&RUN to metoda enzymatyczna służąca do wykrywania oddziaływań białko-DNA w całym genomie. Procedura CUT&RUN rozpoczyna się od związania komórek (lub wyizolowanych jąder komórkowych) z konkanawalną A sprzężoną z kuleczkami magnetycznymi, co umożliwia izolację i manipulację niewielką liczbą komórek w trakcie całej procedury. Wyizolowane komórki poddaje się permeabilizacji przy użyciu łagodnego detergentu, aby ułatwić wprowadzenie przeciwciała skierowanego przeciwko interesującemu białku. Następnie do permeabilizowanych komórek wprowadza się nukleazę mikrokokkową (MNase) połączoną z tagiem białka A lub białka A/G. pA-MNase (lub pAG-MNase) jest rekrutowana do związanego przeciwciała za pomocą tagu białka A lub białka A/G. Po zlokalizowaniu MNase w miejscach docelowych, nukleaza jest krótko aktywowana poprzez dodanie wapnia w celu strawienia DNA wokół białka docelowego. Trawienie przez MNase prowadzi do powstania mononukleosomalnych kompleksów DNA-białko. Następnie wapń jest chelatowany, aby zakończyć reakcję trawienia, a krótkie fragmenty DNA powstałe w wyniku trawienia przez MNase są uwalniane z jąder poprzez krótką inkubację w 37°C, a następnie poddawane oczyszczaniu DNA, przygotowaniu biblioteki i sekwencjonowaniu wysokoprzepustowemu¹. Kliknij tutaj, aby wyświetlić większą wersję tego rysunku.

Rysunek 2: Schematycznego podsumowanie potoku analizy Easy-Shell CUT&RUN. Potok analizy Easy-Shell CUT&RUN został zaprojektowany w trzech głównych sekcjach: (1) kontrola jakości i mapowanie surowych plików odczytów (lewo; kolor fioletowy), (2) normalizacja zmapowanych odczytów oraz liczby odczytów i wyznaczanie piku (środek; kolor zielony) oraz (3) walidacja zmapowanych odczytów i wyznaczonych pików (prawo; kolor różowy). W każdym kroku podano odpowiadający mu numer skryptu powłoki, krótki opis oraz narzędzie programistyczne użyte w tym kroku (w nawiasach). Proste strzałki wskazują bezpośrednie przepływy między krokami. Ten potok analizy CUT&RUN oferuje dwie metody normalizacji odczytów, które mogą spełnić potrzeby użytkowników z odczytami kontrolnymi oraz bez nich, wielowarstwowe procesy walidacji w celu zidentyfikowania odpowiednich powtórzeń do dalszych analiz, a także ukierunkowaną identyfikację pików w celu stworzenia precyzyjnych map ciepła (heatmap) i metaplotów. Potok analizy ten został napisany w formie łatwych w użyciu skryptów powłoki w sposób krok po kroku, aby umożliwić osobom początkującym w bioinformatyce naukę i praktykę podstawowej analizy danych CUT&RUN poprzez czytanie i edytowanie samych skryptów. Kliknij tutaj, aby zobaczyć powiększoną wersję tego rysunku.

Rysunek 3: Porównanie wyników kontroli jakości przed i po przycinaniu jakościowym. Wybrane wyniki raportów kontroli jakości z FastQC pokazują efekt przycinania jakościowego z wykorzystaniem odczytów z SRR8581589 (GSM3609748, CTCF). Przedstawione wyniki obejmują: (A) Wartość jakości dla poszczególnych zasad przed przycinaniem. (B) Ten sam odczyt co w A), po przycinaniu. (C) Rozkład wartości jakości dla wszystkich sekwencji przed przycinaniem. (D) Ten sam odczyt co w C), po przycinaniu. (E) Rozkład GC dla wszystkich sekwencji przed przycinaniem. (F) Ten sam odczyt co w E), po przycinaniu. Minimalna wartość jakości w każdej pozycji w obrębie odczytów sekwencjonowania (A, B) oraz minimalna średnia jakość sekwencji (C, D) ulegają zwiększeniu po przycinaniu jakościowym. Ponadto krok ten może zmniejszyć różnicę między teoretycznym rozkładem liczby zasad GC a rzeczywistą liczbą GC na zasadę w odczytach (E, F) poprzez usunięcie par odczytów o wysokim współczynniku niedopasowania zasad. Kliknij tutaj, aby zobaczyć większą wersję tego rysunku.

Rysunek 4: Analiza rozkładu wielkości wstawek. Histogram wielkości wstawek dla (A) CTCF, (B) H3K27Ac oraz (C) polimerazy RNA II fosforylowanej w serynie 5 (RNAPII-S5P). Histogramy przedstawiają względne różnice w rozkładzie wielkości wstawek między próbkami. Linia przerywana na histogramie reprezentuje skumulowaną frakcję odczytów o wielkości wstawki większej niż lub równej wartości na osi x. n: liczba zgodnych zmapowanych unikalnych odczytów na próbkę po filtracji. FR: fragmenty. Kliknij tutaj, aby zobaczyć powiększoną wersję tego rysunku.

Rysunek 5: Przegląd krajobrazu próbek CUT&RUN. Przedstawiono publicznie dostępne zmapowane odczyty CUT&RUN znormalizowane według przeskalowanego ułamkowego zliczenia (SFRC) bez dodatkowej filtracji (ślady czarne), próbki CUT&RUN przetworzone za pomocą potoku analizy Easy Shells CUTnRUN z normalizacją SFRC (ślady czerwone) oraz „odczyty na milion zmapowanych odczytów w kontroli negatywnej znormalizowane przez spike-in (SRPMC; ślady zielone)” w (A) regionie klastra genów histonów oraz w (B-D) trzech innych regionach z pikami CTCF, H3K27Ac i RNAPII-S5P zidentyfikowanymi przez wszystkie programy do wywoływania pików: MACS2, MACS3 i SEACR. Żółte ramki wyróżniają lokalizację sygnałów spike odfiltrowanych podczas etapu filtracji w potoku analizy Easy Shells CUTnRUN. Kliknij tutaj, aby wyświetlić powiększoną wersję tego rysunku.

Rycina 6: Diagram Venna porównujący piki zidentyfikowane przez różne programy do wywoływania pików (peak callers) oraz różne opcje wywoływania pików. (A) Porównanie pików zidentyfikowanych przez MACS2 i MACS3 z i bez opcji wejściowej IgG podczas wywoływania pików. (B-D) Porównanie pików zidentyfikowanych przez SEACR z i bez opcji wejściowej IgG, z opcjami „stringent” i „relaxed”, oraz z opcją normalizacji przy użyciu plików surowych par odczytów (B), bez opcji normalizacji przy użyciu plików liczb odczytów znormalizowanych metodą SFRC (C) lub plików liczb odczytów znormalizowanych metodą SRPMC (D). (E,F) Porównanie pików zidentyfikowanych przez MACS2, MACS3 z opcją wejściową IgG oraz SEACR z opcją stringent (E) lub relaxed (F). (G,H) Porównanie pików zidentyfikowanych przez SEACR bez opcji wejściowej IgG z opcjami stringent (G) lub relaxed (H). w/ IgG: piki zidentyfikowane z opcją wejściową IgG. w/o IgG: piki zidentyfikowane bez opcji wejściowej IgG. norm: piki zidentyfikowane z opcją normalizacji. non: piki zidentyfikowane bez opcji normalizacji. SFRC: piki zidentyfikowane na podstawie plików liczb odczytów znormalizowanych metodą „scaled fractional count (SFRC)”. SRPMC: piki zidentyfikowane na podstawie plików liczb odczytów znormalizowanych metodą „Spike-in normalized Reads Per Million mapped reads in the negative Control (SRPMC)”. Kliknij tutaj, aby wyświetlić powiększoną wersję tej ryciny.

Rycina 7: Korelacja Pearsona, analiza głównych składowych i diagram Venna w celu walidacji podobieństwa między powtórzeniami. (A-C) Klastrowanie na mapie ciepła z wartościami współczynnika korelacji Pearsona obrazuje stopień podobieństwa między powtórzeniami w obszarach szczytów wyznaczonych przez MACS2 (A), MACS3 (B) oraz SEACR (C). Współczynnik korelacji Pearsona przyjmuje wartości od -1 do 1. Wyższa wartość bezwzględna współczynnika korelacji Pearsona wskazuje na silniejszą korelację między dwiema zmiennymi, a dodatnia wartość wskazuje na korelację dodatnią, w której obie zmienne zmieniają się w tym samym kierunku. W związku z tym próbki o większym podobieństwie wykazują bliższe pokrewieństwo w klastrowaniu na mapie ciepła oraz wyższą wartość współczynnika Pearsona. (D) Analiza głównych składowych (PCA) obrazuje stopień podobieństwa między powtórzeniami i próbkami we wszystkich obszarach szczytów CTCF, H3K27Ac i RNAPII-S5P wyznaczonych przez MACS2 (lewo), MACS3 (środek) i SEACR (prawo). Próbki o większym podobieństwie są rozmieszczone bliżej siebie na wykresie PCA. (E-G) Analiza diagramem Venna w celu porównania szczytów znalezionych w każdym powtórzeniu przez MACS2 (E), MACS3 (F) i SEACR (G). Zaproponowany potok analityczny Easy-Shell CUT&RUN stosuje wszystkie trzy metody w celu zidentyfikowania powtórzeń o wysokim podobieństwie, które mogą być odpowiednie do połączenia wyznaczonych szczytów w dalszych analizach. Proszę kliknąć tutaj, aby wyświetlić powiększoną wersję tej ryciny.

Rysunek 8: Wizualizacja mapy ciepła i metaplotu rozkładu sygnału w pikach. Mapa ciepła i metaploty przedstawiają rozkład wzbogacenia wokół centrów pików wyznaczonych za pomocą różnych programów do wyznaczania pików (peak callers). (A,B) Piki CTCF CUT&RUN wyznaczone z jednej repliki (SRR8581589) przez MACS3 (A) oraz SEACR (B). (C,D) Piki H3K27Ac CUT&RUN wyznaczone z jednej repliki (SRR8581607) przy użyciu MACS3 (C) oraz SEACR (D). (E,F) Piki RNAPII CUT&RUN wyznaczone z jednej repliki (SRR8581589) przez MACS3 (E) oraz SEACR (F). Publicznie dostępne zmapowane pary odczytów ('Public' na Rysunku 8) oraz fragmenty zmapowane przez potok analityczny Easy Shells CUTnRUN ('SFRC' na Rysunku 8) są porównywane po normalizacji metodą 'scaled fractional count (SFRC)'. Piki są wyznaczane przez MACS3 z opcją kontroli IgG ('MACS3 w/ IgG' na Rysunku 8) oraz przez SEACR bez kontroli IgG i bez opcji normalizacji przy użyciu plików liczby odczytów znormalizowanych SFRC w trybie rygorystycznym ('SEACR w/o IgG non SFRC stringent' na Rysunku 8). Przygotowano dwie wersje plików współrzędnych wyznaczonych pików: od początku do końca wyznaczonych pików ('whole' na Rysunku 8) oraz lokalizację bin z najwyższym sygnałem w obrębie wyznaczonych pików (szczyty w pikach wyznaczonych przez MACS3; 'focused' na Rysunku 8). Kliknij tutaj, aby wyświetlić większą wersję tego rysunku.

Tabela 1: Informacje dotyczące plików fastq CUT&RUN w GSE126612. W formie tabeli przedstawiono wszystkie surowe odczyty z plików fastq CUT&RUN zawarte w GSE126612, które zostały wybrane jako przykładowy zestaw danych dla potoku analizy Easy Shells CUTnRUN. Kolumna „File Name” zawiera nazwy plików fastq z surowymi odczytami CUT&RUN, które będą widoczne w katalogu „~/Desktop/GSE126612/fastq” po uruchomieniu skryptu „Script_02_download-fastq.sh”. Kolumna „md5sum” podaje sumy kontrolne MD5 (Message-Digest Algorithm 5) dla przykładowego zestawu danych, które mogą posłużyć do weryfikacji integralności plików po pobraniu zestawu danych za pomocą skryptu „Script_02_download-fastq.sh”. Ostatnia kolumna opisuje cel CUT&RUN dla każdej próbki. Kliknij tutaj, aby pobrać tę tabelę.

Dyskusja

Zdolność do mapowania zajętości białek na chromatynie ma fundamentalne znaczenie dla prowadzenia badań mechanistycznych w dziedzinie biologii chromatyny. W miarę jak laboratoria przyjmują nowe techniki mokrego laboratorium do profilowania chromatyny, możliwość analizowania danych sekwencjonowania z tych eksperymentów w mokrych laboratoriach staje się częstym wąskim gardłem dla naukowców pracujących w mokrych laboratoriach. W związku z tym opisujemy wstępny protokół krok po kroku, aby umożliwić początkującym bioinformatykom pokonanie wąskiego gardła w analizie oraz zainicjowanie analizy i kontroli jakości własnych danych sekwencjonowania CUT&RUN.

Ten protokół analizy CUT&RUN opisuje zastosowanie kilku kroków w celu zapewnienia ilościowego określenia sygnałów w dobrej wierze. Usunięcie odczytów o niskiej jakości i sekwencji adapterów z surowych danych odczytu jest jednym z pierwszych kroków kontroli jakości i jednym z najważniejszych kroków w celu uzyskania dokładnych wyników analizy. W związku z tym ten proces analizy obejmuje łatwe do zastosowania kroki przycinania jakości i adaptera za pomocą programu Trim-galore²⁷. Ze względu na znaczenie tego procesu, ten potok analizy obejmuje etapy porównywania jakości wyników przed (krok 4.3) i po (krok 5.3) procesie przycinania (krok 5.5). Oprócz jakości i przycinania adapterów, ten potok analizy usuwa również niekanoniczne odczyty chromosomów, regiony powtórzeń TA i regiony z czarnej listy, które mogą wprowadzać stronniczość zawartości GC i fałszywie dodatnie skoki/zwane szczyty. Te etapy filtracji stanowią odpowiedni wstęp dla początkujących bioinformatyków, aby zrozumieć krytyczne etapy kontroli jakości w analizie danych CUT&RUN.

Po etapie filtracji ten potok analizy CUT&RUN zapewnia dwie opcje normalizacji: "skalowana liczba odczytów ułamkowych (SFRC)²²" i "Znormalizowane odczyty na milion zmapowane odczyty w kontroli ujemnej (SRPMC)^24,25 w celu utworzenia plików wejściowych do dalszego wywoływania i wizualizacji szczytów. Jeśli oczekuje się, że zestaw danych CUT&RUN ujawni tylko lokalne różnice bez różnic w sygnałach całego genomu między próbkami, skalowana ułamkowa liczba odczytów (ułamek zliczeń pomnożony przez rozmiar referencyjnego krasnala) może wystarczyć do dalszej analizy. Jeśli jednak istnieje możliwość, że między próbkami CUT&RUN wystąpią różnice w sygnale w skali globalnej, użytkownicy mogą wybrać metodę SRPMC, która uwzględnia stosunek odczytów między impulsem a próbką (zarówno eksperymentalne CUT&RUN, jak i próbki kontroli ujemnej) wraz z normalizacją odczytów na milion (RPM) dla odczytów kontroli ujemnej, aby odczyty kontroli ujemnej były porównywalne między różnymi próbkami. Ponieważ SRPMC zapewnia znormalizowane odczyty względem znormalizowanych odczytów kontroli ujemnej, takie podejście minimalizuje sygnał kontroli ujemnej, umożliwiając porównanie zestawów danych utworzonych w różnych partiach i grupach.

Ważnym czynnikiem w wywoływaniu pików próbek CUT&RUN jest eliminacja fałszywie dodatnich pików CUT&RUN podczas analizy in silico , częściowo poprzez włączenie próbek IgG. W szczególności ten potok analizy zapewnia podejścia do wywoływania szczytów dla różnych wywołujących szczyt w celu odrzucenia wyników fałszywie dodatnich CUT&RUN nazywanych szczytami. W przypadku szczytowych wywołujących MACS2/3 nasz potok analizy stosuje próbne odczyty IgG jako próbkę wejściową podczas szczytowego wywołania. W przypadku SEACR w tym procesie analizy zaleca się najpierw niezależne wywołanie pików dla próbek eksperymentalnych i próbek kontroli ujemnej, a następnie usunięcie pików, które nakładają się na próbki eksperymentalne i próbki kontroli ujemnej, ponieważ SEACR może "stracić" większość pików, jeśli zostanie im zapewniona kontrola ujemna podczas wywoływania pików próbek eksperymentalnych. Wyselekcjonowane piki wykazują porównywalne podobieństwo między różnymi szczytowymi rozmówcami i powtórzeniami (Rysunek 5). Ogólnie rzecz biorąc, usunięcie niskiej jakości, niekanonicznych chromosomów, regionu czarnej listy i odczytów powtórzeń TA, przycinanie sekwencji adapterów, normalizacja DNA o skokowym wzroście i właściwa obsługa kontroli negatywnej podczas szczytowych etapów wywoływania zapewnia użytkownikom odpowiednie pliki odczytów, które są odpowiednie do dalszych analiz. Dzięki wysokiej jakości znormalizowanym plikom odczytu i wyselekcjonowanym zwanym szczytami, użytkownicy mogą przystąpić do porównywania podobieństwa między replikacjami oraz tworzyć mapy cieplne i metawykresy z ultraczystym sygnałem tła, aby zweryfikować efektywne wywołanie szczytu.

Wywoływanie szczytów z wysokiej jakości odczytami oznacza rozpoczęcie wyciągania interpretacji biologicznych na podstawie danych CUT&RUN. Protokół ten opisuje uzyskiwanie skoncentrowanych sygnałów szczytowych na mapach cieplnych i metawykresach poprzez wyznaczenie najwyższych lub najbardziej statystycznie istotnych lokalizacji sygnału jako środków pików. Niektóre podejścia do wywołania szczytowego nie wybierają najwyższych sygnałów ani najbardziej istotnych statystycznie sygnałów w ich centralnej lokalizacji. W związku z tym ponowne zdefiniowanie środka każdego piku jako najwyższego sygnału lub najbardziej statystycznie istotnej lokalizacji sygnału służy jako ważny krok do tworzenia wizualnych danych wyjściowych z dobrze skoncentrowanymi sygnałami w środku wykresów. Pliki łóżka oryginalnych nazywanych pików są zachowywane w celu przeprowadzenia adnotacji pików i analizy istotności funkcjonalnej jako kolejne kroki po zakończeniu kroków opisanych w tym protokole.

Chociaż ten potok analizy CUT&RUN zawiera kroki opisujące instalację wymaganych programów, początkujący bioinformatycy mogą napotkać trudności w instalacji narzędzi analitycznych. W związku z tym utworzono powiązaną stronę problemu z Github, aby zapewnić bardziej szczegółowe opisy krok po kroku instalacji programu i ułatwić komunikację w celu wsparcia użytkowników podczas instalacji programu w ich własnym systemie. Kolejne kroki w procesie analizy CUT&RUN, wykraczające poza protokół opisany w tym artykule, obejmują adnotację pików, identyfikację nakładania się różnych typów wywoływanych pików oraz adnotację funkcjonalną dla wywołanych pików. Zakończenie etapów kontroli jakości i wywołanie pików opisane w tym protokole w połączeniu z adnotacją o pikach w dół umożliwi użytkownikom wyciągnięcie biologicznego znaczenia z ich danych CUT&RUN.

Ten potok analizy CUT&RUN został zbudowany w celu zapewnienia ogólnych, wprowadzających wskazówek krok po kroku dotyczących masowej analizy CUT&RUN. Ten potok ma pewne ograniczenia. Po pierwsze, chociaż ten potok analizy próbuje poradzić sobie z efektem zmienności zawartości GC poprzez filtrowanie odczytów w regionach czarnej listy (które obejmują "regiony artefaktów o wysokim sygnale" i "regiony powtórzeń artefaktów" oraz regiony powtórzeń TA), podejście to może nie być wystarczające dla niektórych organizmów, które mogą mieć charakterystyczną zawartość GC w swoim genomie. Dlatego jeśli użytkownicy obawiają się jakichkolwiek uprzedzeń związanych z treścią GC, rozważ dodanie kolejnego kroku w celu poprawienia zmapowanych odczytów. Dla początkujących bioinformatyków "computeGCBias" i "correctGCBias" w Deeptools mogą być opcjami do tego celu. Po drugie, ten potok analizy obsługuje zarówno odczyty o regularnym rozmiarze wkładki (100 bp-1 kb), jak i odczyty o małym rozmiarze wkładki (< 100 pz), które mogą być rzeczywistymi odczytami niektórych białek związanych z chromatyną, w tym samym pliku. Ponieważ ten potok analizy jest napisany w skryptach powłoki, użytkownicy mogą modyfikować "Script_08_bam-to-BEDPE-BED-bedGraph.sh", aby pobierać krótkie odczyty rozmiaru wstawienia oddzielnie od zwykłych odczytów rozmiaru fragmentu podczas etapu generowania zmapowanego pliku łóżka odczytów. Następnie, krótki rozmiar wkładki odczytuje plik łóżka może być znormalizowany niezależnie od zwykłych odczytów mapowanych rozmiarów wkładek, aby zminimalizować efekt zmniejszania. Po trzecie, aby zmniejszyć złożoność procesu analizy, Easy Shells CUTnRUN nie obejmuje etapu próbkowania w dół, aby dopasować głębokość sekwencjonowania próbek CUT&RUN. Użytkownicy mogą jednak zastosować krok próbkowania w dół po przefiltrowaniu plików bam przy użyciu widoku samtools²⁸ lub PositionBasedDownsampleSam (Picard)²⁹.

Wszystkie kroki analizy w tym protokole są napisane w skryptach powłoki, aby umożliwić początkującym bioinformatykom nauczenie się podstaw analizy CUT&RUN poprzez przejrzenie skryptów. Oczekujemy, że użytkownicy będą mogli ćwiczyć analizę bioinformatyczną krok po kroku, uruchamiając każdy skrypt powłoki sekwencyjnie w terminalu. Co więcej, prostota skryptów powłoki dostępnych w tym procesie analizy CUT&RUN pozwala użytkownikom na korygowanie i dostosowywanie tych skryptów w celu zastosowania tego potoku analizy do własnych danych CUT&RUN. Ostatecznie oczekujemy, że ten proces analizy CUT&RUN może zmniejszyć typowe wąskie gardła w procesie analizy danych CUT&RUN, umożliwiając badaczom z mokrych laboratoriów i początkującym bioinformatykom wyciąganie wniosków biologicznych z własnych danych sekwencjonowania CUT&RUN.

Materiały

Lista materiałów użytych w tym artykule

Nazwa	Firma	Numer katalogowy	Komentarze
bedGraphToBigWig	ENCODE	https://hgdownload.soe.ucsc.edu/admin/exe/	Oprogramowanie do kompresji i konwersji readcounts bedGraph do bigWig
bedtools-2.31.1	The Quinlan Lab @ the U. of Utah	https://bedtools.readthedocs.io/en/latest/index.html	Oprogramowanie do przetwarzania plików bam/bed/bedGraph
bowtie2 2.5.4	Uniwersytet	Johnsa Hopkinsahttps://bowtie-bio.sourceforge.net/bowtie2/index.shtml	Oprogramowanie do budowy indeksu muszki i wykonywania wyrównania
CollectInsertSizeMetrics (Picard)	Broad institute	https://github.com/broadinstitute/picard	Oprogramowanie do analizy rozkładu wielkości płytki
Cutadapt	NBIS	https://cutadapt.readthedocs.io/en/stable/index.html	Oprogramowanie do przycinania adapterów
Deeptoolsv3.5.1	Instytut Maxa Plancka	https://deeptools.readthedocs.io/en/develop/index.html	Oprogramowanie do wykonywania analizy korelacji współczynników Pearsona, analizy głównych składowych oraz analizy mapy cieplnej/wykresu średniej
FastQC Wersja 0.12.0	Babraham Bioinformatics	https://github.com/s-andrews/FastQC	Oprogramowanie do sprawdzania jakości pliku fastq
Interwencjav0.6.1	Biologia obliczeniowa i Regulacja genów - grupa Matheliera	https://intervene.readthedocs.io/en/latest/index.html	Oprogramowanie do wykonywania analizy diagramu Venna przy użyciu plików szczytowych
MACSv2.2.9.1	Inicjatywa Chana Zuckerberga	https://github.com/macs3-project/MACS/tree/macs_v2	Oprogramowanie do wywoływania szczytów
MACSv3.0.2	Inicjatywa Chana Zuckerberga	https://github.com/macs3-project/MACS/tree/master	Oprogramowanie do wywoływania szczytów
Samtools-1.21	Wellcome Sanger Institute	https://github.com/samtools/samtools	Oprogramowanie do przetwarzania plików sam/bam
SEACRv1.3	Howard Hughes Medial institute	https://github.com/FredHutch/SEACR	Oprogramowanie do wywoływania szczytów
Toolkit Release 3.1.1	NCBI	https://github.com/ncbi/-tools	Oprogramowanie do pobierania SRR z GEO
Trim_Galore v0.6.10	Babraham Bioinformatics	https://github.com/FelixKrueger/TrimGalore	oprogramowanie do wykonywania wysokiej jakości i szybkiego przycinania