Kwantyfikacja dynamiki cytoszkieletu za pomocą dynamicznej mikroskopii różnicowej

Summary

Dynamiczna mikroskopia różnicowa (DDM) łączy w sobie cechy dynamicznego rozpraszania światła i mikroskopii. W tym artykule przedstawiono proces wykorzystania DDM do scharakteryzowania odtworzonych sieci cytoszkieletów poprzez ilościowe określenie subdyfuzyjnej i zamkniętej w klatce dynamiki cząstek w sieciach wimentynowych oraz ruchu balistycznego aktywnych kompozytów aktynowo-mikrotubulowych napędzanych miozyną.

Abstract

Komórki mogą pełzać, samoleczyć się i dostosowywać swoją sztywność dzięki niezwykle dynamicznemu cytoszkieletowi. W związku z tym odtworzenie sieci biopolimerów cytoszkieletu może prowadzić do powstania wielu aktywnych i elastycznych materiałów. Jednak inżynieria takich materiałów o precyzyjnie dostrojonych właściwościach wymaga zmierzenia, w jaki sposób dynamika zależy od składu sieci i metod syntezy. Kwantyfikacja takiej dynamiki jest wyzwaniem ze względu na różnice w czasie, przestrzeni i przestrzeni formułowania sieci złożonych. Niniejszy protokół opisuje, w jaki sposób technika analizy Fouriera, różnicowa mikroskopia dynamiczna (DDM), może określić ilościowo dynamikę sieci biopolimerowych i jest szczególnie dobrze przystosowana do badań sieci cytoszkieletów. DDM pracuje nad sekwencjami czasowymi obrazów uzyskanych przy użyciu szeregu modalności mikroskopowych, w tym konfokalnego skanowania laserowego, fluorescencji szerokokątnej i obrazowania w jasnym polu. Z takich sekwencji obrazów można wyodrębnić charakterystyczne czasy dekoracji fluktuacji gęstości w rozpiętości wektorów falowych. Opracowano również przyjazny dla użytkownika pakiet Pythona o otwartym kodzie źródłowym do przeprowadzania analizy DDM. Dzięki temu pakietowi można mierzyć dynamikę znakowanych składników cytoszkieletu lub osadzonych cząstek znacznikowych, jak pokazano tutaj na podstawie danych dotyczących pośrednich sieci włókien (wimentyny) i aktywnych sieci aktyna-mikrotubule. Użytkownicy, którzy nie mają wcześniejszego doświadczenia w programowaniu lub przetwarzaniu obrazów, będą mogli wykonać DDM przy użyciu tego pakietu oprogramowania i powiązanej dokumentacji.

Introduction

Cytoszkielet to sieć włókien białkowych, która rozciąga się w poprzek cytoplazmy komórek eukariotycznych, łącząc powierzchnię komórki z jądrem. Ma unikalne właściwości materiałowe, zapewniając mechaniczną ochronę przed dużymi i powtarzającymi się obciążeniami mechanicznymi, a jednocześnie napędzając dynamiczne zmiany kształtu komórki¹. Odtworzone sieci cytoszkieletów mogą dawać początek wielu interesującym dynamicznym zachowaniom, od umieszczania osadzonych cząstek w klatkach po ruch balistyczny napędzany silnikami molekularnymi^{2,3,4,5,6,7,8,9,10,11}. Metody analizy dynamiki takich sieci obejmują śledzenie ruchu wbudowanych mikrosfer znacznikowych^6,7,12,13,14, analiza obrazu w celu śledzenia wielkości klastrów o dużej gęstości białka w czasie⁸, dynamiczne rozpraszanie światła¹⁵, prędkość obrazu cząstek^{4,16,17,18,19}, obliczanie gęstości widmowej mocy obrazów w czasie¹⁹, oraz analiza kymograficzna²⁰. W miarę jak przeprowadzanych jest coraz więcej badań nad odtworzonymi sieciami cytoszkieletów, czy to w celu zrozumienia mechaniki komórkowej, czy materii aktywnej, coraz bardziej potrzebne są solidne, bezstronne i powtarzalne metody charakteryzowania dynamiki. Różnicowa mikroskopia dynamiczna (DDM)^21,22, stosunkowo nowa technika, która została wykorzystana do badania dynamiki cytoszkieletu, jest jedną z takich technik, która skutecznie określa ilościowo dynamikę przy użyciu kilku parametrów zdefiniowanych przez użytkownika. Dzięki opisanemu tutaj pakietowi oprogramowania naukowcy z niewielkim doświadczeniem w programowaniu lub analizie obrazu będą mogli wykorzystać DDM we własnej pracy.

DDM to technika analizy obrazu służąca do wyodrębnienia dynamiki próbki. Podobnie jak śledzenie cząstek lub prędkość obrazu cząstek, DDM wymaga serii czasowej obrazów (często tysięcy obrazów), zwykle rejestrowanych za pomocą mikroskopu. W przeciwieństwie do śledzenia cząstek, poszczególne cechy lub koraliki śledzące nie muszą być zlokalizowane (ani nawet możliwe do zlokalizowania) na obrazie. W przeciwieństwie zarówno do śledzenia cząstek, jak i prędkości obrazu cząstek, za pomocą DDM odzyskuje się dynamikę zespołu przy stosunkowo niewielu parametrach określonych przez użytkownika. Za pomocą DDM obrazy są analizowane w przestrzeni Fouriera w celu określenia czasu zaniku fluktuacji gęstości w zakresie liczb falowych, q, gdzie q = 2πu, a u jest wielkością częstotliwości przestrzennych, . Uzyskuje się informacje podobne do rozpraszania, ale z obrazami w przestrzeni rzeczywistej uzyskanymi pod mikroskopem^21,22,23. W związku z tym można skorzystać z różnych metod mikroskopii generujących kontrast, takich jak fluorescencja szerokiego pola^22,24, fluorescencja konfokalna²⁵, polarized²⁶, dark-field²⁷, lub fluorescencja jasnego arkusza²⁸ mikroskopie. Ponadto obrazy wykorzystywane do analizy DDM mogą być wykorzystywane do śledzenia cząstek lub prędkości obrazu cząstek w celu uzyskania informacji uzupełniających.

Ta kombinacja cech dynamicznego rozpraszania światła i mikroskopii optycznej sprawia, że DDM jest potężną i wszechstronną techniką. Od czasu pierwszego opisu przez Cerbino i Trappe w 2008 roku²¹, gdzie wykazano, że DDM mierzy dyfuzję cząstek koloidalnych o długości fali 73 nm, DDM jest używany do pomiaru przepływających koloidów²⁹, agregacja koloidalna^30,31, lepkosprężystość ciekłych kryształów nematycznych²⁶, dynamika żeli koloidalnych³², pianki gruboziarniste³³, nanocząstki w zamkniętych środowiskach^34,35,36,37, ruchliwość bakterii^38,39,40,41, dyfuzja słabo rozpraszających się klastrów białkowych⁴², fale kapilarne na granicy płynów⁴³ i inne systemy. Osoby poszukujące bardziej kompletnej listy publikacji wykorzystujących DDM mogą zapoznać się z obszernymi pracami przeglądowymi na ten temat^22,23,44,45.

DDM został również wykorzystany do zbadania dynamiki sieci biologicznych. Drechsler i wsp. użyli DDM do pomiaru dynamiki aktyny w żywych oocytach Drosophila⁴⁶. Burla i wsp. określili ilościowo dynamikę cząstek znacznikowych w sieciach kompozytów hialuronanu i hialuronowonowo-kolagenu⁴⁷. Udokumentowano również kilka zastosowań DDM do badania dynamiki cząstek znacznikowych w odtworzonych sieciach cytoszkieletów^9,10, transport cząsteczek DNA w takich sieciach^48,49, a dynamika aktywnych odtworzonych sieci została również udokumentowana^11,50,51. Zaletą DDM w pomiarze dynamiki w takich układach jest to, że poszczególne cząstki lub cząsteczki nie muszą być lokalizowane i śledzone. Tak więc, na przykład, dynamikę cząsteczek DNA w zatłoczonych środowiskach można zmierzyć za pomocą DDM pomimo trudności w śledzeniu tak małych i niesferycznych cząsteczek. Co więcej, dzięki mikroskopii fluorescencyjnej można użyć wielokolorowego znakowania do selektywnego pomiaru dynamiki poszczególnych składników w złożonym kompozycie.

Aby wykonać DDM, wykonywana jest sekwencja obrazów w czasie, I(x,y,t). Dla danego czasu opóźnienia, Δt, znajdowane są wszystkie (lub podzbiór) par obrazów oddzielonych tym czasem opóźnienia. Podniesiona do kwadratu transformata Fouriera różnicy każdej pary,

jest obliczane i uśredniane razem. Wielkość ta, , jest uśredniona promieniowo, pod warunkiem, że dynamika jest izotropowa. W ten sposób powstaje macierz DDM (nazywana również funkcją struktury obrazu), . Proces ten przedstawiono graficznie na rysunku 1. Aby określić dynamikę próbki na podstawie tej macierzy DDM, przyjmuje się, że macierz DDM przyjmuje postać

gdzie A to amplituda, która zależy od szczegółów mikroskopu i struktury próbki, B to tło, które zależy od szumu na obrazach, a f (q, Δt) to pośrednia funkcja rozpraszania (ISF), która zawiera informacje o dynamice^21,22. W prostych przypadkach

gdzie τ to charakterystyczny czas rozpadu lub dekoracji. Taki ISF został wykorzystany w kilku badaniach wykorzystujących DDM w systemach ergodycznych, takich jak rozcieńczone zawiesiny koloidalne^{21,24,27,37,40,52}. Jednak inne formy ISF mogą być wykorzystywane do modelowania różnych typów dynamiki. Na przykład, można użyć ekspansji kumulacyjnej do modelowania ISF dla próbek polidyspersyjnych jako

gdzie μ jest miarą polidyspersyjności^42,53; jeśli fluktuacje gęstości zanikają w dwóch oddzielnych trybach, można użyć ISF takiego jak

^{26, 54, 55, 56, 57};

inne ISFy mogą być używane do pływania mikroorganizmów lub innych aktywnych cząstek^{38,39,40,41,58,59}.

Rysunek 1: Przegląd analizy DDM. Z szeregów czasowych obrazów oblicza się transformację Fouriera różnic obrazów w celu obliczenia macierzy DDM. Macierz DDM można dopasować do modelu w celu określenia skali czasowej wahań gęstości w zakresie wartości q. Kliknij tutaj, aby zobaczyć większą wersję tego rysunku.

Tutaj opisane jest użycie pakietu oprogramowania do analizy DDM opracowanego w Pythonie, PyDDM. Ten pakiet oprogramowania opiera się na pracy wykonanej przez nasze laboratoria badawcze i innych opublikowanych badaniach w ciągu ostatnich kilku lat. Głównymi motywacjami do stworzenia tego pakietu oprogramowania są potrzeba (1) śledzenia i przechowywania metadanych i parametrów używanych w analizie; (2) szczegółową dokumentację ze szczegółowymi przykładami analizy od początku do końca; oraz (3) łatwy sposób na zastosowanie różnych (lub stworzenie nowych) modeli matematycznych do dopasowania danych (np. dodanie modeli ISF, takich jak te niedawno opracowane dla aktywnych filamentów⁶⁰, byłoby proste). Istnieją również inne pakiety oprogramowania do analizy DDM, choć nie wszystkie są dobrze udokumentowane i napisane w języku programowania typu open source. Na przykład istnieje kod C++ z obliczeniami na procesorach graficznych (https://github.com/peterlu/ConDDM)²⁵, kod C++, który wykorzystuje transformacje Fouriera w czasie, aby przyspieszyć obliczenia (https://github.com/giovanni-cerchiari/diffmicro)⁶¹, wersje MATLAB i Python (https://github.com/MathieuLeocmach/DDM)⁴⁰, kod MATLAB (https://sites.engineering.ucsb.edu/~helgeson/ddm.html)²⁷oraz kod MATLAB z kwantyfikacją niepewności (https://github.com/UncertaintyQuantification/DDM-UQ)⁶². Ponieważ ten pakiet PyDDM jest dobrze udokumentowany i zapewnia dużą elastyczność w sposobie obliczania i analizowania macierzy DDM, miejmy nadzieję, że może być przydatny dla badaczy, którzy chcą wdrożyć DDM, niezależnie od ich doświadczenia w programowaniu lub analizie obrazu.

Protokół pokazuje, jak ten pakiet oprogramowania może być wykorzystany do ilościowego określenia dynamiki sieci cytoszkieletów odtworzonych in vitro. Odbywa się to za pomocą dwóch odrębnych zestawów danych obrazowych: (1) obrazów submikronowych cząstek znacznikowych osadzonych w sieci wimentyny wykonanych za pomocą mikroskopii jasnego pola oraz (2) obrazów znakowanych fluorescencyjnie włókien aktyny i mikrotubul w splątanej sieci kompozytowej o aktywności sterowanej miozyną, wykonanych za pomocą laserowej skaningowej mikroskopii konfokalnej. Analizy tych dwóch zestawów danych podkreślają znaczące mocne strony DDM, w tym jego zdolność do analizowania obrazów wykonanych za pomocą różnych modalności obrazowania (np. jasnego pola lub fluorescencji konfokalnej), do wyodrębniania dynamiki z osadzonych znaczników lub z znakowanych włókien oraz do ilościowego określania różnorodności dynamiki (np. subdyfuzyjnej i ograniczonej lub balistycznej).

Protocol

UWAGA: Plik Jupyter Notebook zawierający kod do wykonania każdego kroku w następującym protokole można znaleźć w następującym repozytorium GitHub, https://github.com/rmcgorty/PyDDM/tree/main/Examples. Plik PDF z tym plikiem znajduje się w pliku uzupełniającym 1. Dodatkowo przewodnik po kodzie oraz dokumentację każdej funkcji i klasy można znaleźć na stronie internetowej https://rmcgorty.github.io/PyDDM/.

1. Instalacja oprogramowania

1. Aby postępować zgodnie z przykładowymi plikami analizy DDM, zainstaluj Jupyter Notebook w celu uruchomienia kodu. Zainstaluj inne wymagane typowe pakiety języka Python, w tym NumPy i Matplotlib. Wszystkie te pakiety są dostarczane w pakiecie z dystrybucją Anaconda (zobacz https://www.anaconda.com/products/individual).
2. Zainstaluj pakiet Pythona xarray⁶³. Ten pakiet jest niezbędny do organizowania i przechowywania metadanych i parametrów analizy. Jeśli używasz dystrybucji Anaconda, zainstaluj xarray (wraz z zalecanymi zależnościami) za pomocą polecenia:
   conda install -c conda-forge xarray dask netCDF4 wąskie gardło
3. Zainstaluj pakiet PyYAML za pomocą polecenia:
   conda install -c anaconda yaml
   Ten pakiet jest niezbędny do odczytywania metadanych o obrazach do analizy oraz parametrów ustawionych przez użytkownika do analizy i dopasowania.
4. Zainstaluj pakiet PyDDM, pobierając go z repozytorium GitHub lub używając polecenia git:
   git clone https://github.com/rmcgorty/PyDDM.git

2. Planowanie sesji obrazowania

1. Wybierz optymalną dostępną metodę obrazowania i ustawienia optyczne. Jak wspomniano, DDM może być stosowany z wieloma metodami mikroskopowymi.
2. Aby pomóc w zaplanowaniu odpowiedniej soczewki obiektywu i rozmiaru obrazu do użycia, określ zakres liczb falowych, q, które będą badane na podstawie rozmiaru piksela i całkowitego rozmiaru obrazu. Na podstawie tych obliczeń potwierdź, że wybór powiększenia i pola widzenia jest optymalny dla eksperymentu. W przypadku analizowanych tutaj obrazów użyto obiektywu 60x 1,4 NA i rozmiaru obrazu 256 x 256 pikseli przy rozmiarze piksela 0,83 μm dla aktywnej sieci kompozytowej aktyna-mikrotubula. Do obrazów koralików osadzonych w sieci wimentyny użyto obiektywu 100x 1,4 NA i rozmiaru obrazu 512 x 512 pikseli przy rozmiarze piksela 0,13 μm.
   UWAGA: Minimalne q jest ustawione przez 2π/NΔx, gdzie rozmiar obrazu (zakłada się, że jest kwadratowy) wynosi N × N pikseli przy rozmiarze piksela Δx. Maksimum q to minimum π/Δx i 2π NA/λ, gdzie NA to apertura numeryczna obiektywu do obrazowania, a λ to długość fali światła (w przypadku obrazowania w jasnym polu można zastąpić NA (_obiektyw NA +_kondensor NA)/2).
3. Następnie rozważ zakres skal czasowych do zbadania. Zazwyczaj analiza DDM jest wykonywana na sekwencjach składających się z co najmniej 1000 klatek.
   1. Aby określić odpowiednią liczbę klatek na sekundę, należy wziąć pod uwagę oczekiwany czas potrzebny obiektom w próbie na przemieszczenie się o odległość w kolejności minimalnej skali długości rozpoznawania (odpowiadającej maksymalnej długości q).
   2. Biorąc pod uwagę górną granicę zakresu badanych skal czasowych, należy pamiętać, że zazwyczaj widmo mocy setek różnic obrazów w danym czasie opóźnienia Δt jest uśredniane razem, aby uzyskać wystarczające statystyki do redukcji szumu. W związku z tym należy uzyskać sekwencje obrazów dłuższe niż maksymalna badana skala czasu.
      UWAGA: Jeśli znany jest oczekiwany współczynnik dyfuzji, D lub prędkość v, można oszacować oczekiwane czasy charakterystycznego zaniku za pomocą τ = 1/Dq² lub τ = 1/vq wraz z zakresem q, który został określony na podstawie pola widzenia i rozmiaru piksela. Zakres oczekiwanych wartości τ w dostępnym zakresie q może pomóc w wyborze liczby klatek na sekundę i liczby klatek do uzyskania.

3. Przygotowanie próbki i akwizycja obrazu

UWAGA: Aby uzyskać szczegółowe informacje na temat ustawień przygotowania próbki i obrazowania używanych dla danych przedstawionych w sekcji reprezentatywnych wyników, zobacz poprzednie publikacje autorów^11,51,64 i Plik uzupełniający 2.

1. Opierając się na uwzględnieniu skal czasu i długości do sondowania, uzyskaj sekwencje obrazów, najlepiej ponad 1000 klatek.
   UWAGA: Kod przeanalizuje kwadratowe obrazy lub kwadratowe obszary zainteresowania na obrazie, więc odpowiednio dostosuj rozmiar ramki.
2. Zapisywanie sekwencji obrazów jako trójwymiarowego stosu TIFF w skali szarości. Alternatywnie, format używany przez systemy Nikon Instruments, format ND2, może być odczytywany przez zainstalowany pakiet. Jeśli obrazy są zapisane w innym formacie, użyj programu ImageJ lub innego programu do przetwarzania obrazów, aby przekonwertować obrazy na stos TIFF.
   UWAGA: W przypadku korzystania z plików ND2 musi być zainstalowany pakiet nd2reader z https://github.com/Open-Science-Tools/nd2reader
.

4. Ustawienia parametrów

1. Utwórz kopię pliku parametrów, example_parameter_file.yml znajduje się w repozytorium kodu PyDDM w folderze examples. Otwórz ten plik YAML za pomocą edytora tekstu, takiego jak NotePad++ lub edytora tekstu w JupyterLab. Zobacz Plik uzupełniający 2, aby zapoznać się z przykładowym plikiem parametrów YAML używanym w analizie danych przedstawionych w sekcji wyników reprezentatywnych.
2. W skopiowanym pliku YAML podaj katalog danych i nazwę pliku odpowiadającą sekwencji obrazu, która ma być analizowana. W sekcji metadanych podaj rozmiar w pikselach i liczbę klatek na sekundę.
3. W sekcji Analysis_parameters podaj szczegółowe informacje na temat sposobu obliczania macierzy DDM. Niektóre parametry w tym miejscu są opcjonalne.
   1. Podaj co najmniej wartości parametrów number_lag_times i last_lag_time. Odpowiadają one liczbie różnych czasów opóźnienia, dla których należy obliczyć macierz DDM i najdłuższy czas opóźnienia (w ramkach), który należy zastosować. Dla danych kulek znacznikowych w zastosowanych tu sieciach wimentynowych parametry number_lag_times i last_lag_time wynosiły odpowiednio 60 i 1000. Kod obliczy macierz DDM dla czasów opóźnienia od 1 klatki (lub innego minimalnego czasu opóźnienia, jeśli określono opcjonalny parametr first_lag_time) do last_lag_time z odstępami logarytmicznymi.
      UWAGA: Jeśli uzyskano M-ramki, można by obliczyć macierz DDM dla czasu opóźnienia tak dużego jak M-1. Jednak przy słabych statystykach przy tak dużym opóźnieniu dane mogą być zaszumione. Najdłuższy czas opóźnienia, dla którego należy obliczyć matrycę DDM, będzie zależał od szczegółów danych, ale sugerujemy wypróbowanie około jednej trzeciej całkowitego czasu trwania serii obrazów.
4. W sekcji Fitting_parameters należy podać szczegółowe informacje na temat tego, w jaki sposób należy dopasować macierz DDM lub funkcję rozpraszania pośredniego (ISF). Podaj nazwę modelu pod parametrem modelu. Podaj początkowe przypuszczenie, dolną i górną granicę dla każdego z parametrów dopasowania w wybranym modelu.
   UWAGA: Aby wyświetlić listę możliwych modeli mocowań, uruchom funkcję print_fitting_models. Modele można również znaleźć w dokumentacji online na stronie PyDDM.

5. Obliczanie macierzy DDM

1. Zainicjuj instancję klasy DDM_Analysis. W tym celu podaj metadane i parametry analizy omówione powyżej, przekazując nazwę pliku YAML wraz z pełną ścieżką pliku YAML do DDM_Analysis. Alternatywnie przekaż metadane i parametry jako strukturę danych słownika języka Python.
2. Uruchom funkcję calculate_DDM_matrix, aby obliczyć macierz DDM. To obliczenie może potrwać kilka minut lub dłużej, w zależności od rozmiaru klatki i liczby opóźnień. Zobacz Rysunek 2, aby zapoznać się z typowymi czasami działania.
3. Sprawdź zwrócone dane, które będą znajdować się w strukturze danych z pakietu xarray znanego jako zestaw danych. Ta struktura danych jest przechowywana pod atrybutem ddm_dataset.
   UWAGA: W tej strukturze danych będzie przechowywana nie tylko macierz DDM, ale także powiązane zmienne i metadane. Zostanie on również zapisany na dysku w formacie Network Common Data Form (netCDF).
4. Sprawdź wykresy i rysunki, które zostaną wygenerowane i wyświetlone. Rysunki te są również zapisywane jako plik PDF w katalogu danych.
   1. Zobacz, że jeden z wygenerowanych wykresów pokazuje uśredniony przez zespół kwadratowy moduł obrazów przekształconych przez Fouriera, jako funkcja q. Domyślnie kod używa tego do oszacowania parametru tła B. Oszacuj tło z zakładając, że w granicy dużego q zbliży się do B/2, gdzie B jest tłem.
   2. Jeśli nie osiąga plateau przy dużym q, użyj innej metody szacowania B. Aby to osiągnąć, ustaw parametr background_method w pliku YAML lub jako opcjonalny argument słowa kluczowego do funkcji calculate_DDM_matrix. Więcej szczegółów na temat metod szacowania B przedstawiono w sekcji dotyczącej wyników reprezentatywnych.

Rysunek 2: Czas obliczeń dla obliczenia macierzy DDM. W (A) i (B) czasie obliczania macierzy DDM pokazane jest . Dane używane we wszystkich przypadkach to film składający się z 5000 klatek o rozmiarze obrazu 512 x 512 pikseli. Macierz DDM została obliczona dla 30 czasów opóźnienia, logarytmicznie rozmieszczonych między 1 klatką (0,01 s) a 1000 ramek (10 s). Kod został uruchomiony na komputerze stacjonarnym z procesorem Intel i7-10700 2,90 GHz i 32 GB pamięci RAM. W (A) pokazano efekt różnicowania liczby różnic obrazów używanych do obliczania macierzy DDM dla każdego czasu opóźnienia. W tym celu obrazy są łączone do binów, aby uzyskać rozmiar obrazu 256 x 256. Dla każdego czasu opóźnienia Δt obrazy oddzielone tym Δt są odejmowane, a otrzymana macierz jest przekształcana Fourierem. Dla danego Δt można użyć wszystkich par obrazów oddzielonych tym Δt (pokazanych na niebiesko), można użyć tylko nienakładających się par obrazów (np. klatki 1 i 10, 10 i 19 itd.; pokazane na brązowo) lub można użyć 300 par obrazów lub mniej dla każdego Δt. W (B) pokazany jest wpływ zmiany rozmiaru obrazu na czas obliczeń. Obrazy zostały podzielone na grupy przez zgrupowanie 2 x 2, 4 x 4 lub 8 x 8 pikseli, co dało rozmiary obrazów odpowiednio 256 x 256, 128 x 128 lub 64 x 64. Dla każdego z nich używa się około 300 par obrazów do obliczania macierzy DDM dla każdego Δt. (C) Z macierzy DDM można wyodrębnić funkcję rozpraszania pośredniego (ISF). Jest to pokazane dla trzech przypadków w (A). Niebieskie punkty danych (bez przesunięcia) odpowiadają ISF, gdy dla każdego Δt stosuje się maksymalną liczbę par obrazów; brązowe punkty danych (z przesunięciem 0,1) odpowiadają ISF, gdy dla każdego Δt stosuje się nienakładające się na siebie pary obrazów; a różowe punkty danych (z przesunięciem 0,2) odpowiadają ISF, gdy dla każdego Δt używanych jest co najwyżej 300 par obrazów. ISF znaleziony przy użyciu nienakładających się par obrazów wykazuje hałaśliwość przy długości Δt. W takim przypadku używa się kilku par obrazów przy długim Δt (np. dla Δt 1000 klatek używane są tylko 4 pary obrazów). (D) Poprzez dopasowanie ISF do funkcji wykładniczej określa się charakterystyczny czas zaniku, τ, dla każdej liczby falowej q. W kolorze różowym wyniki są wyświetlane po skażeniu oryginalnych obrazów w formacie 2 x 2, co daje rozmiar obrazu 256 x 256. W kolorze szarym wyniki są wyświetlane po kategoryzacji przez 8 x 8, co daje rozmiar obrazu 64 x 64. Poprzez kategorowanie danych informacje o dynamice przy wyższych liczbach falowych są tracone, ale obliczanie macierzy DDM dla obrazów 64 x 64 jest około 16 razy szybsze niż dla obrazów 256 x 256. Kliknij tutaj, aby zobaczyć większą wersję tego rysunku.

6. Dopasowanie matrycy DDM lub ISF

1. Zainicjuj instancję klasy DDM_Fit. W tym celu należy przekazać do DDM_Fit nazwę pliku YAML zawierającego metadane obrazu i parametry dopasowania.
2. Zdecyduj, który model macierzy DDM lub ISF ma być używany do dopasowywania danych. Wyświetl listę dostępnych modeli, wykonując funkcję print_fitting_models. Określ model, który ma być używany w pliku parametrów YAML lub przy użyciu funkcji reload_fit_model_by_name.
3. Ustaw początkowe domysły i granice dla każdego parametru w wybranym modelu w podanym pliku parametrów YAML. Aby zmienić początkowe odgadnięcie dla dowolnego parametru, użyj funkcji set_parameter_initial_guess. Ustaw granice parametrów za pomocą funkcji set_parameter_bounds. Na przykład, jak widać w pliku uzupełniającym 2, dla danych kulek śledzących w sieci wimentyny, początkowe przypuszczenie czasu rozpadu wynosiło 1 s, a granice tego parametru wynosiły 0,01 s i 2000 s.
4. Wykonaj dopasowanie za pomocą dopasowania funkcyjnego. Przypisz zmienną do danych wyjściowych tej funkcji, aby łatwo uzyskać dostęp do wyników.
   UWAGA: Ta funkcja może przyjmować wiele opcjonalnych argumentów. Zapoznaj się z dokumentacją kodu i podanymi przykładami, aby uzyskać listę takich argumentów i kiedy należy rozważyć ustawienie ich na wartości inne niż domyślne.

7. Interpretacja wyników dopasowania

1. Generuj wykresy do kontroli pasowań i zależności q parametrów dopasowania za pomocą funkcji fit_report.
   UWAGA: Ta funkcja wygeneruje serię wykresów, które zostaną również zapisane jako plik PDF. Opcjonalne argumenty tej funkcji mogą być używane do modyfikowania utworzonych wykresów.
2. Wśród wygenerowanych wykresów znajdzie się rysunek z wykresami cząstkowymi 2 x 2 przedstawiający macierz DDM lub ISF (w zależności od wybranego modelu dopasowania) przy czterech wartościach q (określonych jako opcjonalny argument do fit_report), wraz z obliczoną macierzą DDM lub ISF przy użyciu parametrów modelu i najlepszego dopasowania. Aby wykreślić macierz DDM lub ISF wraz z najlepszym dopasowaniem w interaktywny sposób, użyj Browse_DDM_Fits klasy, jak pokazano w podanych przykładach, gdy jest używane środowisko Jupyter Notebook.
3. Na podstawie wykresu charakterystycznego czasu zaniku τ w funkcji liczby falowej q określ, czy dynamika podąża za ruchem dyfuzyjnym, subdyfuzyjnym, balistycznym czy innym rodzajem ruchu. Można to zrobić, szukając relacji prawa potęgowego między τ a q.
   UWAGA: Na wykresie logarytmicznym τ w funkcji q wygenerowanym przez funkcję fit_report zostaną pokazane trzy linie, odpowiadające dopasowaniom prawa potęgowego w określonym zakresie wartości q. Ciągła linia odpowiada dopasowaniu τ vs. q do prawa potęgowego, τ = 1 / Kq^β, gdzie K i β są parametrami dowolnymi. Linia przerywana w kolorze pomarańczowym odpowiada dopasowaniu do prostej dyfuzji, τ = 1 / Dq², gdzie D jest współczynnikiem dyfuzji. Kropkowana linia w kolorze niebieskim odpowiada dopasowaniu do τ = 1 / vq, gdzie v jest prędkością.

8. Zapisywanie wyników

1. Wyniki dopasowania zostaną zapisane w zestawie danych xarray. Użyj funkcji xarray to_netcdf lub wbudowanego modułu pickle Pythona, aby zapisać tę strukturę danych na dysku. Użyj funkcji xarray open_dataset, aby załadować te pliki netCDF.
2. Użyj funkcji save_fit_results_to_excel, aby zapisać wyniki dopasowania wraz z danymi w pliku arkusza.

Results

Tutaj pokazujemy przykłady analizy wykonanej za pomocą PyDDM z dwóch różnych zestawów eksperymentów. W jednym zestawie eksperymentów submikronowe kulki znacznikowe zostały osadzone w sieciach składających się z pośredniego białka filamentu wimentyny i zobrazowane za pomocą soczewki obiektywowej 100x w trybie jasnego pola przy 100 klatkach / s (Rysunek 3A). Wimentyna ulega ekspresji w komórkach mezenchymalnych i jest kluczowym wyznacznikiem właściwości mechanicznych cytoplazmy⁶⁵ oraz stabilności mechanicznej jądra w komórkach wykonujących migrację ograniczoną^66,67. Do tej pory zrekonstruowane sieci wimentyny były badane głównie za pomocą reologii makroskopowej^64,68,69, podczas gdy dynamice poświęcono stosunkowo mało uwagi^13,70,71. Dodatkowe informacje na temat tych eksperymentów można znaleźć w pliku uzupełniającym 2. W drugim zestawie eksperymentów przygotowano aktywne sieci cytoszkieletów z aktyną, mikrotubulami i miozyną. Spektralnie wyraźne znaczniki fluorescencyjne pozwoliły na obrazowanie filamentów aktyny i mikrotubul za pomocą dwukolorowego laserowego skaningowego mikroskopu konfokalnego przy użyciu 60-krotnej soczewki obiektywowej przy 2,78 klatki/s (Rysunek 3B, C). Aktyna i włókna mikrotubul są ważnymi czynnikami napędzającymi dynamiczne zmiany kształtu komórek, a ich działania są koordynowane przez interakcje mechaniczne i biochemiczne⁷². Dodatkowe informacje na temat tych eksperymentów można znaleźć w¹¹. Poszczególne klatki z sekwencji obrazów wykonanych w tych eksperymentach są pokazane w Rysunek 3.

Rysunek 3: Analizowane obrazy z szeregów czasowych. (A) Obraz jasnego pola z koralikami 0,6 μm w sieci wimentyny. (B,C) Obraz (B) mikrotubul i (C) aktyny w aktywnym kompozycie aktyna-mikrotubula wykonany obiektywem 60x na laserowym skaningowym mikroskopie konfokalnym, przy użyciu światła wzbudzającego 561 nm do obrazowania mikrotubul i światła wzbudzającego 488 nm do obrazowania aktyny. Kliknij tutaj, aby zobaczyć większą wersję tego rysunku.

Dla obrazów koralików śledzących w sieciach vimentin, nagrano filmy o długości 5000 klatek i rozmiarze 512 x 512 pikseli przy 100 klatkach na sekundę. Na ich podstawie obliczono macierz DDM przy 60 logarytmicznie rozmieszczonych czasach opóźnienia od 1 do 1000 klatek, czyli od 0,01 s do 10 s. Aby oszacować tło, B, średnia kwadratowych obrazów przekształconych metodą Fouriera, , obliczono i ustawiono na wartość ^55,73. Przyjęto założenie, że w ciągu największych 10% wartości q wielkość ta jest równa B/2 i że B jest niezależne od q. Jest to domyślna metoda szacowania B w pakiecie, ale inne metody są możliwe, ustawiając parametr background_method na inną wartość.

Z parametrami A(q) i B wyznaczonymi na podstawie , można wyodrębnić funkcję rozpraszania pośredniego (ISF) z macierzy DDM. Przykładowe ISF są pokazane w Rysunek 4. W Rysunek 4A, ISF z obrazów kulek o średnicy 0,6 μm osadzonych w sieci o stężeniu wimentyny 19 μM jest pokazany. W Rysunek 4B, pokazany jest ISF dla tego samego typu koralików w sieci o stężeniu wimentyny 34 μM. Co ciekawe, w żadnym z tych przypadków ISF nie spadł do zera. Przy dużych opóźnieniach ISF powinien zbliżać się do zera dla systemów ergodycznych. Oznacza to, że w takich systemach wahania gęstości powinny być całkowicie dekoracyjne w przypadku dużych czasów opóźnienia. Fakt, że ISF nie rozpadł się tutaj do zera, mógł wynikać z niedokładnych oszacowań A(q) i B, które posłużyły do znalezienia ISF na podstawie obliczonej macierzy DDM. Warto zauważyć, że zastosowana tutaj metoda może zawyżać B w niektórych scenariuszach⁶². Jednak bardziej prawdopodobne jest, że dynamika kulek śledzących jest naprawdę nieergodyczna, ponieważ koraliki mają rozmiar porównywalny do rozmiaru siatki sieci i dlatego mogą zostać zamknięte w klatce. Inne dane potwierdziły stwierdzenie nieergodyczności. Mianowicie, rozmiar kulki 0,6 μm był większy niż obliczona średnia wartość dla rozmiarów oczek wynosząca 0,4 μm dla stężenia 19 μM i 0,3 μm dla stężenia 34 μM. Dodatkowo, wyniki śledzenia pojedynczych cząstek tych kulek znacznikowych, które są pokazane później, również wykazały ograniczony ruch.

Rysunek 4: Pośrednie funkcje rozpraszania przy kilku liczbach falowych dla sieci vimentin. ISF jest wykreślany jako funkcja czasu opóźnienia dla wartości q od około 1 do 9 ^μm-1. (A) ISF z obrazów kulek 0,6 μm w sieci wimentyny o stężeniu wimentyny 19 μM. (B) ISF z obrazów kulek 0,6 μm w sieci wimentyny o stężeniu wimentyny 34 μM. Długie plateau czasowe ISF przy wartości znacznie powyżej zera wskazuje na nieergodyczność. Kliknij tutaj, aby zobaczyć większą wersję tego rysunku.

Biorąc pod uwagę, że dynamika jest prawdopodobnie nieergodyczna, ISF są dopasowane do postaci , gdzie C jest czynnikiem nieergodyczności ³². Ta forma ISF była używana w poprzednich badaniach dynamiki nieergodycznej, takich jak żele koloidalne^32,74 lub cząstki znacznikowe w sieciach aktyna-mikrotubule ¹⁰. Kropkowane czarne linie w Rysunek 4 pokazują dopasowania wraz z danymi. Na podstawie tych dopasowań można teraz spojrzeć na zależność q czasu rozpadu τ i parametru nieergodyczności C.

Rysunek 5: Czas zaniku w funkcji liczby falowej dla sieci vimentin. Od pasowań do ISF, czas zaniku τ jest określany dla zakresu wartości q. Dla jasności nie pokazujemy wartości τ dla każdego q, ale tylko zbiór logarytmicznie rozstawiony. W kolorze niebieskim (jasnobrązowym) znajdują się dane z obrazów kulek 0,6 μm w sieciach wimentyny o stężeniu wimentyny 19 μM (34 μM). Słupki błędów reprezentują odchylenia standardowe w τ dla wielu filmów (cztery filmy dla danych z siecią 19 μM [niebieski] i pięć filmów dla danych z siecią 34 μM [tan]). Czerwone linie przerywane wyznaczają szacunkowe granice naszej rozdzielczości czasowej i przestrzennej, zgodnie z opisem w wynikach. Ciągła linia pokazuje skalowanie, co wskazywałoby na ruch dyfuzyjny. Żaden zestaw danych nie jest zgodny z tym skalowaniem. Koraliki w sieci 19 μM wykazują raczej ruch subdyfuzyjny ( z >β > 2), a koraliki w sieci 34 μM wykazują ruch ograniczony lub w klatce. Kliknij tutaj, aby zobaczyć większą wersję tego rysunku.

Czasy rozpadu wykazały dużą dozę niepewności, zarówno na niskim jak i wysokim ekstremum q, co widać na Rysunek 5. Słupki błędów na tym wykresie pokazują odchylenie standardowe między czterema filmami analizowanymi dla przypadku niższego stężenia wimentyny lub pięcioma filmami analizowanymi dla wyższego stężenia. Aby zrozumieć źródło dużej niepewności w tych skrajnościach, należy wziąć pod uwagę zarówno rozdzielczość czasową, jak i przestrzenną. Przybliżone granice rozdzielczości są pokazane trzema czerwonymi liniami przerywanymi. Dwie poziome linie odpowiadają minimalnemu i maksymalnemu czasowi opóźnienia sondowanego. Biorąc pod uwagę liczbę klatek na sekundę wynoszącą 100 klatek/s i maksymalny czas opóźnienia odpowiadający 1000 klatek (20% całkowitego czasu trwania wideo), dokładność została utracona podczas pomiaru dynamiki zachodzącej szybciej niż 0,01 s lub wolniej niż 10 s. Przy niższych wartościach q dopasowane wartości dla τ były większe niż 10 s. W związku z tym należy spodziewać się dużych niepewności w przypadku czasów zaniku, które są większe niż maksymalny czas opóźnienia. Na górnym końcu zakresu q czas zaniku zbliżał się do minimalnego czasu opóźnienia wynoszącego 0,01 s, ale utrzymywał się powyżej tego poziomu. Zamiast być ograniczonym przez rozdzielczość czasową, przy tych wyższych wartościach q, rozdzielczość przestrzenna może być czynnikiem ograniczającym. Biorąc pod uwagę rozmiar piksela wynoszący 0,13 μm, największa wartość q wynosiła około 24 ^μm-1. Jednak rozdzielczość ograniczona dyfrakcją niekoniecznie pozwala na dokładne pomiary dynamiki przy tak wysokich częstotliwościach przestrzennych. Przybliżenie rozdzielczości optycznej jako prowadzi do górnej granicy liczby falowej około 16 ^μm-1, biorąc pod uwagę aperturę numeryczną soczewki, NA, wynoszącą 1,4 i długość fali światła, . Jest to oznaczone pionową czerwoną przerywaną linią w Rysunek 5. Rzeczywiście, dane były zaszumione przy dużych wartościach q. Jeszcze przed tą przybliżoną górną granicą q zaobserwowano zwiększoną niepewność τ, co mogło wynikać z przeszacowania q_max. Gorsza rozdzielczość optyczna niż przewidywano może wynikać z tego, że do obrazowania poza szkiełkiem nakrywkowym do próbki wodnej użyto soczewki olejowej lub dlatego, że soczewka kondensatora była nieidealnie wyrównana.

Dla koralików 0,6 μm osadzonych w mniej skoncentrowanej sieci (19 μM wimentyna), można zaobserwować na wykresie logarytmicznym czasu zaniku w funkcji liczby falowej, że czas zaniku zmniejszał się wraz z liczbą falową w sposób zgodny z prawem potęgowym (Rysunek 5). Wydaje się jednak, że nie jest to zgodne z tym, czego można by oczekiwać dla normalnego ruchu dyfuzyjnego, gdzie . Raczej τ zmniejszało się bardziej stromo wraz ze wzrostem q. Wskazuje to na ruch subdyfuzyjny, który często występuje w przypadku koralików w zatłoczonych środowiskach, takich jak te. Dopasowanie τ(q) w zakresie od 1,4 ^μm-1 do 12,3 ^μm-1 do prawa mocy w postaci τ = 1/Kq^β daje parametry transportu K = 0,0953^{μm β} / s i β = 2,2. Dla tych, którzy są bardziej przyzwyczajeni do myślenia o normalnej dyfuzji i subdyfuzji w kategoriach średniego kwadratu przemieszczenia (MSD) cząstek znacznikowych w funkcji czasu opóźnienia (tj. MSD = K' Δt^α), pomocne jest uznanie, że wykładnik skalowania subdyfuzyjnego w równaniu MSD, α, jest równoważny α = 2 / β. Innymi słowy, wartość β = 2,2 jest zgodna z wykładnikiem skalowania subdyfuzyjnego w równaniu MSD wynoszącym α = 0,9. Można ustawić PyDDM tak, aby pasował τ(q) do tego zakresu wartości q, określając indeksy tablicy q za pomocą parametru Good_q_range w pliku YAML lub przekazując opcjonalny argument forced_qs do funkcji generate_fit_report. Zakres q od 1,4 ^μm-1 do 12,3 ^μm-1 odpowiadałby dla danych tutaj indeksom tablicy q od 15 do 130.

Dla koralików 0,6 μm w bardziej skoncentrowanej sieci (34 μM), czas rozpadu wykazywał niewielką zależność od q. Jest to prawdopodobnie spowodowane nieergodycznością koralików w sieci o mniejszym rozmiarze oczek. Aby zbadać nieergodyczność w tym systemie, parametr nieergodyczności, C, powinien zostać wykreślony jako funkcja q, jak w Rysunek 6. Dla kulek 0,6 μm w sieci 19 μM wimentyny, C ≈ 0,2 z niewielką zależnością od q (nie pokazano). Jednak dla sieci z 34 μM wimentyną i dla sieci o jeszcze wyższym stężeniu 49 μM wimentyny, log C był proporcjonalny do q², jak pokazano na Rysunek 6. Ta zależność między C i q jest oczekiwana dla ruchu ograniczonego. W przypadku koralików uwięzionych w kieszeniach sieci, oczekuje się, że MSD ustabilizuje się w wystarczająco długim czasie opóźnienia (tj. , gdzie to MSD, a δ² to maksymalny MSD). Ponieważ ISF może być wyrażony w kategoriach MSD jako , a nieergodyczny ISF przechodzi do C z dużym opóźnieniem (tj. ), relacja jest uzyskiwany^32,75. Dlatego można użyć C(q), aby znaleźć δ², a to daje δ² = 0,017 μm² i 0,0032 μm² odpowiednio dla sieci wimentyny 34 i 49 μM (co odpowiada δ = 0,13 μm i 0,057 μm).

Rysunek 6: Parametr nonergodicity vs. liczba falowa dla sieci vimentin. Od pasowań do ISF, parametr nieergodyczności C jest określany dla zakresu wartości q. W kolorze jasnobrązowym (czerwonym) znajdują się dane z obrazów kulek 0,6 μm w sieciach wimentyny o stężeniu wimentyny 34 μM (49 μM). Słupki błędów reprezentują odchylenia standardowe w τ dla wielu filmów (pięć filmów dla danych z siecią 34 μM [tan] i cztery filmy dla danych z siecią 49 μM [czerwony]). Oś y ma skalowanie logarytmiczne. Obserwuje się zależność q od C, która następuje , co pozwala na wyodrębnienie maksymalnego średniego kwadratu przemieszczenia, δ². Pasuje do są oznaczone liniami ciągłymi. Kliknij tutaj, aby zobaczyć większą wersję tego rysunku.

Można użyć innych metod, aby wyodrębnić δ wielkości zamknięcia z danych, jak również wykładnik subdyfuzyjny znaleziony podczas badania τ(q) dla koralików w sieci wimentyny 19 μM. Po pierwsze, można użyć metody opisanej przez Bayles et al.⁷⁶ oraz Edera et al.⁷⁷, aby wyodrębnić MSD z macierzy DDM. Warto zauważyć, że ta metoda nie wymaga dopasowania matrycy DDM. Wystarczy obliczyć macierz DDM, D(q, Δt) i (z którego można wyznaczyć A(q) i B). Następnie, aby znaleźć MSD, należy użyć relacji . Należy zauważyć, że ta metoda znajdowania MSD zakłada, że rozkład przemieszczeń cząstek jest gaussowski, chociaż wcześniejsze prace wykazały, że w niektórych przypadkach zaburzenia układu mięśniowo-szkieletowego pochodzące z DDM zgadzają się z zaburzeniami układu mięśniowo-szkieletowego ze śledzenia cząstek, nawet jeśli przemieszczenia są niegaussowskie⁷³. W przypadku tego systemu, zgodnie z oczekiwaniami⁷⁸, występuje niegaussowska rozkład dużych przemieszczeń, jak widać na rysunku S1. W pakiecie PyDDM powinna zostać wykonana funkcja extract_MSD, która zwraca . Po drugie, można użyć śledzenia pojedynczych cząstek, aby znaleźć MSD. Chociaż DDM może być używany do analizy obrazów, w których albo wysoka gęstość cząstek, albo ograniczona rozdzielczość optyczna uniemożliwia dokładną lokalizację cząstek, w przypadku obrazów kulek 0,6 μm w sieciach wimentynowych byliśmy w stanie zlokalizować i śledzić koraliki za pomocą oprogramowania trackpy (https://github.com/soft-matter/trackpy)⁷⁹. Ten pakiet oprogramowania do śledzenia cząstek wykorzystuje algorytmy opisane przez Crockera i Griera⁸⁰.

Rysunek 7: Średnie kwadratowe przemieszczenie w funkcji czasu opóźnienia dla sieci vimentin. Zaburzenia układu mięśniowo-szkieletowego oznaczono dwiema metodami. Po pierwsze, MSD obliczono z macierzy DDM (pokazanej za pomocą symboli stałych). Następnie określono MSD za pomocą śledzenia pojedynczych cząstek (SPT) w celu znalezienia trajektorii cząstek (symbole otwarte). Słupki błędów są wyznaczane w taki sam sposób, jak opisano w dwóch poprzednich legendach rysunków. (A) Zaburzenia układu mięśniowo-szkieletowego dla kulek 0,6 μm w sieci wimentyny 19 μM wskazują na ruch subdyfuzyjny, z dobrą zgodnością między tymi dwiema metodami znajdowania zaburzeń układu mięśniowo-szkieletowego. (B) Zaburzenia układu mięśniowo-szkieletowego dla kulek 0,6 μm w sieci wimentyny 49 μM wskazują na ruch w klatce, z dobrą zgodnością między dwiema metodami znajdowania MSD i z maksymalną MSD stwierdzoną na podstawie parametru nieergodyczności. Kliknij tutaj, aby zobaczyć większą wersję tego rysunku.

Zaburzenia układu mięśniowo-szkieletowego w porównaniu z czasem opóźnienia dla koralików 0,6 μm w sieci 19 μM i w sieci 49 μM są pokazane w Rysunek 7. W obu przypadkach zaburzenia układu mięśniowo-szkieletowego wyznaczone na podstawie DDM dobrze zgadzały się z zaburzeniami układu mięśniowo-szkieletowego uzyskanymi za pomocą śledzenia pojedynczych cząstek (SPT). Ponadto, dla mniej skoncentrowanej sieci, wykładnik skalowania subdyfuzyjnego (α w ) wynosił około 0,9. Jest to zgodne ze skalowaniem τ(q) znalezione przez dopasowanie ISF do wyznaczenia τ(q) (czyli 2/2,2 = 0,9). W przypadku bardziej skoncentrowanej sieci MSD stabilizuje się przy dłuższych czasach opóźnienia. Maksymalne MSD znalezione przez analizę zależności q parametru nonergodicity (pokazanego na Rysunek 7B z poziomą linią na δ² = 0,0032 μm²) było w przybliżeniu tą samą wartością, do której MSD zarówno z SPT, jak i DDM wydawały się zmierzać. Istnieje rozbieżność między najdłuższymi zaburzeniami układu mięśniowo-szkieletowego o czasie opóźnienia określonymi na podstawie DDM i SPT w Rysunek 7A. Chociaż może to wynikać z ograniczonej liczby długich trajektorii czasu opóźnienia, może być również tak, że dalsza optymalizacja zakresu wartości q, dla których macierz DDM jest używana do oszacowania dla każdego czasu opóźnienia (tak jak zrobili Bayles et al.⁷⁶ i Edera et al.⁷⁷) poprawiłoby nasze wyniki, I taka optymalizacja będzie przedmiotem przyszłych prac.

Te eksperymenty, w których rejestrowano sekwencje obrazów kulek znacznikowych osadzonych w sieci włókien pośrednich wimentyny, pozwoliły na niezależne analizy: DDM (przy użyciu pakietu opisanego tutaj) i SPT (przy użyciu trackpy). Obie analizy mogą ujawnić stopień subdyfuzji i długość uwięzienia, co pozwala na użycie dwóch niezależnych technik analizy obrazu w celu uzyskania uzupełniających się wskaźników. Istnieją dodatkowe ilości, które można porównać z SPT i DDM. Na przykład niejednorodność dynamiki próbki może objawić się jako niegaussowska w rozkładzie przemieszczeń cząstek (tj. rozkładzie van Hove'a) wyznaczonym na podstawie SPT, a także w ISF wyznaczonym na podstawie DDM, który pasuje do rozciągniętego wykładniczego ^34,35. Rysunek S1 przedstawia rozkład van Hove'a dla cząstek 0,6 μm w sieciach wimentyny i omawia wykładnik rozciągania uzyskany po dopasowaniu ISF - metryk używanych w tandemie w poprzednich badaniach w celu zademonstrowania heterogenicznej dynamiki cząstek w systemach biomimetycznych9^,10,47 lub innych zatłoczonych środowisk ³⁴. Innym przykładem jest ISF obliczony na podstawie trajektorii cząstek zmierzonych za pomocą SPT i porównanych z ISF nabytymi przez DDM. Podczas gdy średnie kwadratowe przemieszczenia i rozkłady przemieszczeń są metrykami najczęściej pobieranymi z analizy SPT, można również obliczyć ISF na podstawie trajektorii cząstek, , używając (patrz rysunek S2). Ten ISF może być porównany z ISF generowanymi przez DDM i użyty do ujawnienia dynamiki, która nie jest widoczna w klasie MSD⁵⁹.

Podczas gdy pozyskiwanie obrazów cząstek znacznikowych w sieci może pozwolić na korzystanie z uzupełniających metod analizy SPT i DDM, ważne jest, aby zauważyć, że przewagą DDM nad SPT jest to, że nie wymaga obrazów koralików (lub innych cech), które można łatwo zlokalizować i śledzić. Aby zademonstrować ten punkt, zwracamy następnie uwagę na analizę aktywnych sieci filamentów aktynowych i mikrotubul, w których fluorescencyjne znakowanie aktyny i tubuliny pozwala na obrazowanie obu typów filamentów, różniących się od siebie za pomocą różnych fluoroforów, za pomocą wielokolorowego laserowo-skaningowego mikroskopu konfokalnego.

Obrazy zostały uzyskane za pomocą laserowo-skaningowego mikroskopu konfokalnego sieci aktyna-mikrotubula z aktywnością napędzaną przez miozynę (miozyna mięśni szkieletowych królika II; Cytoszkielet #MY02). Szczegóły eksperymentów i wyników zostały wcześniej opisane ¹¹, a reprezentatywne wyniki pokazane tutaj pochodzą z analizy dwóch filmów dostarczonych w materiałach uzupełniających (filmy S1 i S4) dla¹¹. Obie sekwencje obrazów zostały zarejestrowane z prędkością 2,78 klatki/s przez 1000 klatek.

Aby przeanalizować te obrazy, matryca DDM została obliczona dla 50 czasów opóźnienia w zakresie od 0,4 s do 252 s (1 klatka do 700 klatek). Macierz DDM została następnie dopasowana do modelu , przy czym pośrednią funkcją rozpraszania jest . Istnieją zatem cztery parametry dopasowania: A, τ, s i B. Wyniki tych dopasowań są pokazane w Rysunek 8. Zaobserwowano, że macierz DDM dla określonej wartości q miała plateau w krótkich czasach opóźnienia, wzrastała wraz z czasem opóźnienia, a następnie ustabilizowała się (lub wykazywała oznaki początku plateau) w dużych czasach opóźnienia. Macierz DDM dla niższych wartości q nie osiągnęła plateau przy długich czasach opóźnienia. Należy zatem spodziewać się słabej dokładności w pomiarze czasu zaniku dla tych niskich dynamik q (duża skala długości).

Charakterystyczne czasy zaniku, τ, od dopasowań do macierzy DDM są pokazane w Rysunek 9. Wyniki przedstawiono dla aktywnej sieci złożonej aktyna-mikrotubula (podobnie jak w filmie S1¹¹) oraz dla aktywnej sieci aktynowej (podobnie jak w filmie S4¹¹). Obie sieci zostały przygotowane z tymi samymi stężeniami aktyny i miozyny, ale sieć zawierająca tylko aktynę została utworzona bez tubuliny, jak opisano w¹¹. Dla tych dwóch typów aktywnych sieci obserwowana zależność prawa potęgowego wynosiła . To skalowanie wskazuje na ruch balistyczny i na to, że skurcz i przepływ napędzany miozyną dominuje nad ruchem termicznym włókien. Od τ = (vq)^-1 można znaleźć prędkość charakterystyczną v wynoszącą około 10 nm/s dla aktywnej sieci aktyna-mikrotubula i 75 nm/s dla aktywnej sieci aktyny. Wartości te są zgodne z analizą prędkości obrazu cząstek w tych samych filmach pokazanych w¹¹. Skalowanie skalowanie nie utrzymywało się na niższych wartościach q dla aktywnej sieci kompozytowej aktyna-mikrotubula. Jest to prawdopodobne, ponieważ rzeczywiste czasy rozpadu dla tej złożonej sieci aktynowo-mikrotubulowej przy niższych wartościach q są dłuższe niż maksymalny czas opóźnienia obliczonej matrycy DDM. Maksymalny czas opóźnienia jest oznaczony poziomą czerwoną linią w Rysunek 9, a czasy zaniku odbiegają od oczekiwanego skalując się w pobliżu tych dłuższych czasów.

Rysunek 8: Macierz DDM vs. czas opóźnienia dla aktywnej sieci kompozytowej aktyna-mikrotubula. Macierz DDM dla kilku wartości q jest wykreślana jako funkcja czasu opóźnienia z filmu sieci złożonej złożonej z 2,9 μM monomerów aktyny, 2,9 μM dimerów tubuliny i 0,24 μM miozyny. Dane te pokazują analizę tylko kanału mikrotubul w wielokolorowych szeregach czasowych obrazów. Kliknij tutaj, aby zobaczyć większą wersję tego rysunku.

Rysunek 9: Czas zaniku w funkcji liczby falowej dla aktywnych sieci aktyna-mikrotubule. Z dopasowania macierzy DDM uzyskuje się czas zaniku, τ, w funkcji liczby falowej q. Wykreślone to τ vs q dla obrazów aktywnej sieci aktyna-mikrotubula (analizująca tylko kanał mikrotubul) w kolorze brązowym i dla obrazów aktywnej sieci aktyny w kolorze zielonym. Obie sieci mają takie same stężenia aktyny i miozyny (odpowiednio 2,9 μM i 0,24 μM); kompozyt aktyna-mikrotubula ma 2,9 μM dimerów tubuliny. Czasy rozpadu dla aktywnej sieci aktynowej są znacznie krótsze niż czasy rozpadu dla aktywnej sieci aktyna-mikrotubula, co wskazuje na szybszy ruch aktywnej sieci aktyny. W obu przypadkach dynamika jest balistyczna, ponieważ dane są zgodne z trendem trend. Wstawka: wykres ISF w funkcji czasu opóźnienia przeskalowanego przez liczbę falową (Δt × q) pokazuje załamanie ISF w zakresie wartości q. Wskazuje to również na ruch balistyczny. ISF pokazane w tej wstawce pochodzą z aktywnej sieci aktyny. Kliknij tutaj, aby zobaczyć większą wersję tego rysunku.

Dla tych danych aktywnych sieci, zdecydowaliśmy się dopasować matrycę DDM, . Kontrastuje to z tym, co zrobiono dla danych o koralikach w sieci wimentyny, gdzie A (q) i B zostały oszacowane bez żadnego dopasowania do izolacji ISF, f (q, Δt). W tym przypadku dla aktywnych danych sieciowych A i B pozostawiono jako parametry dopasowania, ponieważ metody użyte do oszacowania B nie dały dobrych dopasowań. Domyślną metodą szacowania B jest obliczenie i założenie, że przy dużym q przechodzi to do B/2. Jednak metoda ta przeszacowała B dla tych danych, co było widoczne w fakcie, że przy obliczaniu ISF z B oszacowanych w ten sposób (nie pokazano), ISF były większe niż 1 we wczesnych czasach opóźnienia (podczas gdy powinny przejść od maksimum 1 do zera lub jakiegoś parametru nieergodyczności wraz ze wzrostem czasu opóźnienia). Można wybrać inne metody szacowania B za pomocą parametru background_method. Jedną z tych innych metod jest oszacowanie B jako minimum macierzy DDM we wczesnych czasach opóźnienia (ustawionych na background_method=1). Podobna metoda została użyta przez Baylesa i wsp.⁷⁶, chociaż nie założyli, że B jest stałe z q. Inną opcją jest oszacowanie B jako średniej wartości ze wszystkich czasów opóźnienia macierzy DDM przy maksymalnym q (ustawionym na background_method=2). Te różne metody szacowania tła, jak również wyniki pozwalające na to, aby B było dowolnym parametrem dopasowania, są pokazane w Rysunek 10. Z tych wykresów widać, że amplituda A nie osiągnęła zera przy największych badanych wartościach q, ponieważ nie ustabilizowała się na dużym q (Rysunek 10B), a ponieważ D(q_max, Δt) przeszła z niższego plateau czasu opóźnienia do pewnego wyższego plateau czasu opóźnienia (tj. przy _{q max}, było niezerowe A; Rysunek 10D). Dlatego ani oszacowanie B jako ani jako byłoby odpowiednie. Należy sprawdzić vs. q i D(q_max, Δt) vs. Δt przed podjęciem decyzji o tym, jak (lub czy) oszacować B.

Rysunek 10: Tło a liczba falowa dla aktywnych sieci aktyna-mikrotubule. Z dopasowania macierzy DDM można znaleźć tło, B, jako funkcję liczby falowej, q. Pokazano B w porównaniu z q dla obrazów aktywnej sieci aktyna-mikrotubula (analizująca tylko kanał mikrotubul) wyznaczonej na podstawie tych dopasowań za pomocą fioletowych symboli. Trzy ciągłe linie w (A) pokazują szacunki tła znalezionego bez żadnego dopasowania. Górna, najciemniejsza linia w (A) pokazuje szacowane tło za pomocą , co może być odpowiednie, jeśli plateaus do stałej wartości na dużym q. Z (B) zauważ, że nie osiągnął jeszcze stałej wartości przy największym badanym q. Dlatego użycie tej metody przeszacowuje tło. Dolna linia w (A) pokazuje szacowane tło za pomocą . Jeśli macierz DDM wykazuje plateau czasu opóźnienia o niskim poziomie, jak pokazano w (C) z czerwoną linią, wówczas ta metoda może być odpowiednia do oszacowania tła. Środkowa, najjaśniejsza linia w (A) pokazuje szacowane tło od . Metoda ta może być właściwa, jeżeli przy _{q max} amplituda A osiągnęła zero. Z (D) widać, że amplituda jest różna od zera i dlatego ta metoda przeszacowuje tło. Kliknij tutaj, aby zobaczyć większą wersję tego rysunku.

Rysunek uzupełniający S1: Rozkłady prawdopodobieństwa przemieszczeń cząstek. Rozkłady prawdopodobieństwa przemieszczeń cząstek wykazują niegaussowskość dla stężeń wimentyny 34 μM i 49 μM. Śledzenie pojedynczych cząstek o średnicy 0,6 μm przeprowadzono w sieciach wimentyny o różnych stężeniach. Różne czasy opóźnienia są pokazane w rozkładach przemieszczeń dla trzech warunków. (A) Rozkład przemieszczeń cząstek w sieci wimentyny 19 μM jest zgodny z funkcją Gaussa. Szerokość gaussa zwiększa się wraz ze wzrostem czasu opóźnienia. (B) Rozkład przemieszczeń cząstek w sieci wimentyny 34 μM wykazuje większą niegaussowskość, zwłaszcza przy dużych przemieszczeniach, niż w przypadku 19 μM. (C) Rozkład przemieszczeń cząstek w sieci wimentyny 49 μM również wykazuje niegaussowskość. Co więcej, szerokości rozkładów nie zwiększają się wraz z czasem opóźnienia tak znacząco, jak w próbkach o niższych stężeniach wimentyny, co wskazuje na ograniczony ruch. Niegaussowskie rozkłady van Hove'a (obserwowane dla wszystkich próbek wimentyny, ale najbardziej widoczne w wyższych stężeniach) są związane z niejednorodną dynamiką, co często obserwuje się w transporcie cząstek w zatłoczonych i ograniczonych środowiskach. Innym wskaźnikiem transportu heterogenicznego, który jest określany na podstawie analizy DDM, jest wykładnik rozciągający używany do dopasowania pośredniej funkcji rozpraszania (parametr s w równaniu dla ISF użyty tutaj: + ). Średnie wykładniki rozciągania w zakresie q od 0,4 ^μm-1 do 9,4 ^μm-1 wynoszą, od najwyższego stężenia wimentyny do najniższego, 0,53 ± 0,07, 0,64 ± 0,02 i 0,86 ± 0,04 (średnia ± odchylenia standardowego). Kliknij tutaj, aby pobrać ten plik.

Rysunek uzupełniający S2: Pośrednie funkcje rozpraszania z DDM i SPT. Pokazano pośrednie funkcje rozpraszania (ISF) dla pięciu różnych liczb falowych. ISF w funkcji czasu opóźnienia znalezionego przez DDM jest wykreślany za pomocą okrągłych znaczników, a ISF jest obliczany na podstawie trajektorii pojedynczych cząstek z otwartymi kwadratami. Przerywane czarne linie pokazują dopasowania do ISF nabytych przez DDM. ISF jest obliczany na podstawie trajektorii pojedynczych cząstek, , używając . W (A) ISF pokazano dla cząstek 0,6 μm w sieciach 19 μM wimentyny. W (B) ISF pokazano dla cząstek 0,6 μm w sieciach wimentyny 34 μM. Rozbieżności w ISF stwierdzone przez DDM i SPT są prawdopodobnie spowodowane ograniczoną liczbą długich trajektorii opóźnień. Kliknij tutaj, aby pobrać ten plik.

Plik uzupełniający 1: Protokół do korzystania z DDM. Prezentowane są dane wejściowe i wyjściowe kroków pokazanych w protokole. Kliknij tutaj, aby pobrać ten plik.

Plik uzupełniający 2: Szczegóły przygotowania próbki i przykładowe pliki parametrów dla sieci vimentin. Przedstawiono szczegółowe kroki przygotowania próbki i akwizycji obrazu w sieciach wimentyny. Dodatkowo udostępniono przykładowy plik parametrów do analizy danych przedstawionych w sekcji reprezentatywnych wyników dotyczących sieci wimentynowych. Kliknij tutaj, aby pobrać ten plik.

Discussion

Opisany tutaj pakiet oprogramowania wykorzystuje DDM do analizy fluktuacji gęstości obserwowanych w obrazach uzyskanych za pomocą mikroskopu optycznego. W pierwszej kolejności przedstawiono reprezentatywne wyniki z danych dotyczących cząstek znacznikowych osadzonych w sieciach wimentyny. Analiza takich danych może być wykorzystana do scharakteryzowania rozmiaru siatki i sztywności sieci, podobnie jak śledzenie pojedynczych cząstek było wykorzystywane w wielu wcześniejszych badaniach sieci cytoszkieletów ^6,12,13. Zaletą stosowania DDM w porównaniu ze śledzeniem pojedynczych cząstek jest to, że DDM nie wymaga lokalizacji cząstek. Dlatego nawet na obrazach, na których gęstość cząstek jest zbyt wysoka lub cząstki są zbyt małe, aby można je było zlokalizować i śledzić, DDM może nadal określać dynamikę. Śledzenie pojedynczych cząstek byłoby korzystne w przypadku kontroli zmienności między cząstkami. W przypadku DDM można znaleźć uśrednioną dynamikę zespołu, podczas gdy przy śledzeniu pojedynczej cząstki można obliczyć zarówno MSD pojedynczej cząstki, jak i MSD uśrednione zespołowo. Jednak DDM może być używany do badania dynamiki heterogenicznej poprzez analizowanie wielu obszarów zainteresowania w dużym polu widzenia.

Następnie przedstawiono reprezentatywne wyniki z danych włókien znakowanych fluorescencyjnie w sieci aktywnej złożonej z dwóch różnie znakowanych typów włókien cytoszkieletu¹¹. Na podstawie tych danych scharakteryzowano ruch balistyczny bez konieczności lokalizowania jakichkolwiek cech na obrazie. Ponieważ DDM ekstrahuje uśrednioną dynamikę zespołu przy niewielkiej liczbie danych wejściowych użytkownika, ułatwia to porównywanie serii obrazów uzyskanych w różnych warunkach (np. porównywanie próbek o różnych proporcjach aktyny do mikrotubul lub próbek o różnych stężeniach miozyny, jak to miało miejsce w⁵⁰). Dodatkowo, korzystając z obrazowania fluorescencyjnego, możemy badać dynamikę różnych elementów sieci za pomocą wielokolorowego znakowania. Dokonano tego w ^11,50, gdzie dynamika aktyny i mikrotubul była oddzielnie analizowana w aktywnej sieci kompozytowej aktyna-mikrotubula przy użyciu obrazowania wielokolorowego. W sekcji wyników reprezentatywnych pokazano tylko wyniki z kanału mikrotubul, ale w poprzedniej pracy porównaliśmy dynamikę włókien mikrotubuli i aktyny¹¹.

Zauważamy, że te reprezentatywne wyniki pokazują albo pasywną subdyfuzję, albo aktywny ruch balistyczny. Co ważne, DDM może być wykorzystywany do analizy układów, w których występuje przecięcie w rodzaju dynamiki w pośrednich skalach czasowych lub długościowych. Jako przykłady, Kurzthaler i in. wykorzystali DDM z systemem aktywnych koloidów Janusa do zbadania aktywnego ruchu kierowanego w krótkich skalach czasowych i randomizacji orientacji w dłuższych skalach czasowych⁵⁹; Giavazzi i wsp. użyli DDM z gruboziarnistą pianką i znaleźli skrzyżowanie w dynamice odpowiadającej skali długości pęcherzyka³³; a Cho i in. użyli DDM z żelami koloidalnymi i znaleźli trzy rozróżnialne reżimy w różnych skalach długości, rozciągających się od klastrów fraktalnych do całej sieci³².

Dane zawarte w sekcji wyników reprezentatywnych uzyskano za pomocą mikroskopii jasnego pola i laserowo-skaningowej mikroskopii konfokalnej. Jednak, jak wspomniano wcześniej, DDM może być stosowany z wieloma modalnościami obrazowania. W przypadku każdej metody obrazowania użytkownicy powinni wziąć pod uwagę ustawienia optyczne, takie jak stopień przekroju optycznego lub głębia ostrości. Wysoki stopień przekroju optycznego może zmniejszyć sygnał od nieostrych obiektów, ale nie będzie w stanie dokładnie zmierzyć dynamiki w skalach czasowych większych niż skala czasu, w której obiekty poruszają się poza głębią ostrości^25,28. Bardziej szczegółowe omówienie tego, jak głębia ostrości zależna od q wpływa na analizę DDM, można znaleźć w²². W przypadku obrazowania w jasnym polu użytkownicy mogą być również zmuszeni do rozważenia grubości próbki. Podczas gdy w przypadku próbek słabo rozpraszających, grubsze próbki mogą dostarczyć więcej sygnału⁴², próbki mętne mogą wymagać modyfikacji analizy w celu uwzględnienia rozpraszania wielokrotnego⁸¹. Wreszcie, w przypadku metod obrazowania, które nie są niezmiennicze w przestrzeni liniowej (to znaczy, gdy intensywność zarejestrowana przez kamerę obiektu zależy od tego, gdzie ten obiekt znajduje się w płaszczyźnie próbki x-y), może być konieczne uwzględnienie liniowej wariancji przestrzeni, jak wykazano w przypadku DDM²⁷ w ciemnym polu.

Dla tych, którzy rozpoczynają pracę z DDM, pragniemy podkreślić znaczenie rozważenia rozdzielczości przestrzennej i czasowej. Podczas sprawdzania wyznaczonych czasów zaniku w funkcji liczby falowej ważne jest, aby zaznaczyć granice swojej rozdzielczości (tj. maksymalne i minimalne czasy opóźnienia oraz maksymalną liczbę falową, jak pokazano na rysunku 5). Przed zebraniem danych należy dokładnie przemyśleć te ograniczenia, aby można było wybrać optymalny obiektyw, rozmiar obrazu, liczbę klatek na sekundę i czas trwania filmu. Inną ważną kwestią jest sposób oszacowania parametru tła B. W literaturze przedmiotu stosowano wiele metod szacowania tła, a skutki przeszacowania lub niedoszacowania B opisano we wcześniejszych publikacjach^62,77. Jak pokazano na rysunku 10, PyDDM umożliwia użytkownikom implementację różnych metod szacowania B i sugerujemy, aby nowi użytkownicy wypróbowali te metody i ocenili, które są odpowiednie do użycia.

Mocną stroną tego pakietu jest jego dokładna dokumentacja i przewodniki z przykładowymi danymi, przechowywanie i organizacja metadanych w celu śledzenia sposobu przeprowadzania analiz oraz elastyczność w sposobie analizy macierzy DDM (różne modele dopasowania, wiele metod szacowania parametru tła B, możliwość znalezienia MSD). Istnieje jednak wiele aspektów tego kodu, które można ulepszyć. Obecnie kod nie został zoptymalizowany pod kątem dużej szybkości obliczeń. Metody przyspieszające obliczenia zostały zgłoszone^61,62 i zostaną zaimplementowane w przyszłych wersjach. Ponadto planujemy wdrożyć ostatnio zgłoszone metody, aby lepiej oszacować niepewności i wykorzystać symulacje do kierowania użytkowników do odpowiedniego modelu ISF⁶². Mamy nadzieję, że w przypadku innych ulepszeń użytkownicy skontaktują się z nami z sugestiami.

Materials

List of materials used in this article

Name	Company	Catalog Number	Comments
Kamera CMOS, Orca-Flash 4.0	Hamatsu
F-127, Pluronic	Sigma, Aldrich
Jupyter, Notebook
Nanodrop	, Thermo Fisher
Nikon Ti-Eclipse, mikroskop	Nikon
PLL-PEG-bio	SuSos AG, Dü bendorf, Szwajcaria
Koraliki polistyrenowe	Sigma Aldrich
Minikaseta do dializy białkowej	Thermo Fisher
PyDDM	University of San Diego	N/A Oprogramowanie	open source dostępne w https://github.com/rmcgorty/PyDDM