Artykuł metodologiczny

Ilościowa zależność struktura-aktywność, przewidywanie aktywności i dynamika molekularna nienukleotydowych inhibitorów odwrotnej transkryptazy

2.4K wyświetleń

⸱

DOI:

10.3791/67457

⸱

9 maja 2025

 ,  ,  , 
1Department of Chemical Sciences, Faculty of Science, University of Johannesburg

W tym artykule

Podsumowanie

W tym badaniu wykorzystano strategie in-silico, aby zidentyfikować enumerowaną etrawirynę jako obiecujący środek terapeutyczny dla HIV. Nasze odkrycia dotyczące interakcji molekularnych i dynamiki wspierają racjonalne projektowanie nowych NNRTI jako możliwych alternatyw leczenia HIV.

Streszczenie

Rosnąca częstość występowania lekooporności HIV-1 stanowi wyzwanie dla skuteczności skojarzonej terapii antyretrowirusowej, szczególnie w Afryce Południowej. Rozwój oporności na nienukleotydowe inhibitory odwrotnej transkryptazy (NNRTI) zagraża długoterminowemu powodzeniu terapii przeciwretrowirusowej. Szacuje się, że w 2019 r. oporność na środki przeciwdrobnoustrojowe była bezpośrednio przyczyną 1,27 mln zgonów na całym świecie. W badaniu tym zastosowano podejście in-silico w celu zbadania leków NNTRI i ich pochodnych. Stosowane techniki obejmowały obliczenia teorii funkcjonału gęstości, dokowanie molekularne, wyliczanie, ilościową analizę zależności struktura-aktywność (QSAR), symulację dynamiki molekularnej (MDS) oraz mechanikę molekularną z uogólnionymi metodami Borna i powierzchni. Analiza koncentrowała się na różnych pochodnych pirymidyny i sześciu lekach NNRTI, badając ich interakcje z białkiem HIV-1 (kod PDB 1HQU).

Model QSAR został opracowany w celu przewidywania aktywności biologicznej sześciu badanych NNRTI. Wykorzystując 94 pochodne pirymidyny, model QSAR osiągnął R2 NA poziomie 0,822 i Q2 na poziomie 0,815, co wskazuje na wysoki poziom dokładności predykcyjnej.

MDS przeprowadzono w celu oceny stabilności różnych ligandów i ich nowo opracowanych alternatyw, upewniając się, że pozostają one związane z miejscem aktywnym białka przez 200-nanosekundowy okres symulacji. Etrawiryna wykazywała wahania średniego odchylenia kwadratowego (RMSD) wynoszące około 4,5 A, podczas gdy jej wyliczone pochodne wykazywały wahania RMSD wynoszące 3,5 A. Dzięki dokowaniu molekularnemu, MDS i obliczeniom energii swobodnej, enumeracyjna etrawiryna wykazała najlepszą wydajność, z wartością aktywności 7,373 i wynikiem dokowania -10,517 kcal/mol. Co więcej, obliczona energia swobodna wiązania dla enumerowanej etrawiryny wynosiła -89,684 kcal/mol, przewyższając inne badane ligandy. Znacząca poprawa sugeruje, że zmodyfikowana etrawiryna ma obiecujący potencjał jako nowy środek w terapii antyretrowirusowej.

Uzyskana niższa wartość RMSD, zwiększone interakcje aminokwasów i najwyższa wolna energia wiązania wskazują, że wymieniona etrawiryna może służyć jako realna alternatywa w leczeniu HIV/AIDS.

Wprowadzenie

Pomimo znacznych postępów w leczeniu, poważnym globalnym zagrożeniem dla zdrowia, jakim jest zespół nabytego niedoboru odporności (AIDS), spowodowany przez ludzki wirus niedoboru odporności-1 (HIV-1), pozostaje stałym zagrożeniem1. Leki przeciwretrowirusowe znane jako inhibitory odwrotnej transkryptazy, RTI, są stosowane w leczeniu zakażenia wirusem HIV. Odwrotna transkryptaza, wirusowy enzym polimerazy kwasu dezoksyrybonukleinowego (DNA) niezbędny do replikacji retrowirusa, jest hamowana przez inhibitory odwrotnej transkryptazy (RTI). Nukleozydowe inhibitory odwrotnej transkryptazy (NRTI) i nienukleotydowe inhibitory odwrotnej transkryptazy (NNRTI) są głównymi RTIs2.

Wysoce aktywna terapia antyretrowirusowa (HAART), połączenie różnych leków przeciwwirusowych, stała się standardową terapią HIV, skutecznie kontrolując rozprzestrzenianie się AIDS i przekształcając tę niegdyś śmiertelną chorobę w możliwą do opanowania chorobę przewlekłą3. NNRTI dla HIV-1 są obecnie istotną częścią schematu HAART4. Oporność na środki przeciwdrobnoustrojowe (AMR) jest jednym z najpoważniejszych globalnych zagrożeń dla zdrowia, z szacunkową liczbą 1,27 miliona ofiar śmiertelnych5. W związku z tym istnieje pilna potrzeba rozwiązania problemu oporności na środki przeciwdrobnoustrojowe i zidentyfikowania nowych środków przeciwdrobnoustrojowych. Według Światowej Organizacji Zdrowia (WHO) oporność na środki przeciwdrobnoustrojowe osiągnęła alarmujący poziom w wielu częściach świata.

Stanowi to poważne zagrożenie dla osiągnięcia celów zrównoważonego rozwoju, ponieważ podważa bezpieczeństwo żywnościowe, wzrost gospodarczy i bezpieczeństwo zdrowotne, a także przyczynia się do nierówności społecznych i ekonomicznych6. Częstość występowania wirusów lekoopornych rośnie, nawet w lekach antyretrowirusowych przeznaczonych do zwalczania wirusa HIV. Według najnowszych statystyk na koniec 2022 roku prawie 30 milionów ludzi na całym świecie przyjmowało terapię antyretrowirusową7.

WHO przeprowadziło 30 ankiet i odkryło, że w 21 z tych badań ponad 10% osób, które rozpoczęły terapię antyretrowirusową pierwszego rzutu, miało oporność na Newirapinę lub Efawirenz8. Co więcej, osoby, które wcześniej miały kontakt z lekami antyretrowirusowymi, są do trzech razy bardziej narażone na oporność na NNRTI niż osoby bez ekspozycji9. Badania wykazały, że znaczna liczba niemowląt w wieku poniżej 18 miesięcy, u których niedawno zdiagnozowano HIV, wykazywała wysokie wskaźniki szczepów lekoopornych10,11. Szokujące jest to, że prawie połowa z nich miała szczep oporny na - (NNRTI) jeszcze przed rozpoczęciem leczenia. Odkrycia te podkreślają potrzebę przyspieszenia badań w celu opracowania innowacyjnej terapii HIV.

Rozwój szczepów opornych na leki i niepożądane skutki uboczne wynikające z długotrwałego stosowania nieuchronnie stanowiły wyzwanie dla klinicznego zastosowania NNRTIs12. Rozpoczęcie skojarzonej terapii przeciwretrowirusowej (cART) wydłuża oczekiwaną długość życia osób żyjących z HIV pomimo potencjalnych działań niepożądanych. Te działania niepożądane obejmują ryzyko rozwoju chorób niezakaźnych, w tym lipodystrofii, hiperlipidemii, zmniejszonej gęstości mineralnej kości, podwyższonego poziomu glukozy we krwi prowadzącego do cukrzycy typu 2, nadciśnienia, zwiększonego ryzyka udaru mózgu i problemów związanych z otyłością13. Wczesna diagnoza i szybki dostęp do odpowiedniej opieki medycznej w początkowej fazie zakażenia wirusem HIV oferują znaczne korzyści zarówno z perspektywy klinicznej, jak i publicznej. Szybkie rozpoczęcie ART i profilaktyka przeciwko zakażeniom oportunistycznym prowadzi do znacznego zmniejszenia liczby chorób i śmiertelności związanej z HIV.

Użycie ART może również przyczynić się do zmniejszenia możliwości dalszego przenoszenia wirusa HIV poprzez zmniejszenie poziomu krążącego kwasu rybonukleinowego HIV. Co więcej, zajęcie się leczeniem innych chorób przenoszonych drogą płciową i koinfekcji może również zmniejszyć prawdopodobieństwo dalszego zakażenia wirusem HIV14. NNRTI należą do opcjonalnych klas leczenia inhibitorów, które wchodzą w interakcje z RT poprzez wiązanie się z regionem allosterycznym lub miejscem zakażenia wirusem HIV. Ten rodzaj hamowania jest powszechnie znany jako inhibitor niekonkurencyjny, ponieważ NNRTI nie wiąże się w miejscu aktywnym substratu, ale raczej na zewnątrz. Efekt ten zmienia konformację miejsca wiązania substratu, utrudniając wiązanie standardowego podłoża i prowadząc do przedwczesnego zakończenia łańcucha. Ze względu na ich zmniejszoną toksyczność w porównaniu z NRTI, prostą strukturę, lepszą biodostępność w porównaniu z inhibitorami proteazy i doskonałą selektywność, NNRTI stały się najbardziej atrakcyjnymi inhibitorami HIV15. W związku z tym synteza i projektowanie nowych NNRTI mają kluczowe znaczenie z perspektywy farmakokinetycznej16,17.

NNRTI są niezbędne w leczeniu HIV/AIDS ze względu na ich silną skuteczność i niską toksyczność. Jednak wczesne NNRTI, takie jak newirapina, delawirdyna i efawirenz, napotykają oporność na mutacje wirusa w miejscu wiązania NNRTI18. Leki te, należące do rodziny diarylopirymidyny (DAPY), wykazują silną aktywność przeciwko różnym szczepom NNRTIs, w tym tym opornym na wczesne NNRTIs19, prawdopodobnie ze względu na ich elastyczność molekularną i wyższą barierę odporności na HIV-120,21. Pomimo sukcesu, wysoki wskaźnik mutacji w HIV-1 RT i brak wewnętrznej aktywności korektorskiej doprowadziły do powstania nowych profili oporności u pacjentów stosujących etrawirynę i rylpiwirynę22,23. Te mutacje wirusowe różnią się od siebie, podobnie jak te związane z wczesnymi lekami NNRTI24. Ponad 50 strukturalnie zróżnicowanych klas związków zostało zidentyfikowanych jako NNRTI. Warto zauważyć, że sześć NNRTI uzyskało zgodę na leczenie HIV-1. Do tych zatwierdzonych leków należą: newirapina (NVP), delawirdyna (DLV), efawirenz (EFV), etrawiryna (ETR), rylpiwiryna (RPV) i dorodyryna (DOR). Rysunek 1 pokazuje struktury chemiczne tych sześciu zatwierdzonych leków NNRTI25.

W niedawnym badaniu26, obliczenia DFT i dokowanie molekularne sugerują, że lowastatyna i symwastatyna wykazują potencjał jako środki przeciw koronawirusowi. Wirtualne badania przesiewowe pozwoliły zidentyfikować pięć nowych zatwierdzonych przez FDA cząsteczek kandydujących podobnych do rusztowania efawirenzu, o doskonałym powinowactwie wiązania w aktywnej kieszeni głównej proteazy COVID-19.

Soltan et al.27 przeprowadził podobne badanie nad identyfikacją nowych cząsteczek w celu poprawy zdolności wiązania z HIV RT, stosując strategię opartą na fragmentach z lekami zatwierdzonymi przez FDA do projektowania pochodnych chemicznych. W szczególności wykorzystali struktury delawirdyny, efawirenzu, etrawiryny i rylpiwiryny jako podstawowe rusztowania. Podobieństwo tych pochodnych do leków oceniono za pomocą Swiss-ADME, a następnie zadokowano je w powiązanych strukturach krystalicznych. Badanie zakończyło się wyborem związków wykazujących większe powinowactwo wiązania w porównaniu z ich macierzystymi rusztowaniami, w szczególności zwracając uwagę na wyraźniejszą poprawę w przypadku pochodnych zaprojektowanych z NNRTI drugiej generacji, etrawiryny i rylpiwiryny. Na przykład pochodne RPV01 i RPV15 wykazały znaczną poprawę wartości energii zadokowanej w porównaniu z rylpiwiryną, co wskazuje na potencjalną użyteczność w zwalczaniu zarówno dzikich, jak i zmutowanych form HIV RT.

Murugesan i współpracownicy28 przeprowadzili badania, proponując różne podejścia chemii medycznej w celu zwiększenia skuteczności i zminimalizowania oporności w leczeniu HIV. W badaniu wykorzystano hybrydyzację molekularną, wymianę bioizosteryczną i wysokoprzepustowe badania przesiewowe. W ramach swoich badań udało im się zidentyfikować nowe rusztowania NNRTI o wysokiej sile działania zarówno przeciwko dzikim, jak i lekoopornym szczepom wirusa HIV. Opracowali pochodne DAPY, które wykazują doskonałą selektywność i niską toksyczność, a niektóre związki wykazują skuteczne hamowanie przy stężeniach nanomolowych.

Ostatnio, używanie narzędzi obliczeniowych wraz z badaniami in-silico zyskało popularność dzięki praktycznej analizie właściwości chemicznych leków kwantowych29. W tym badaniu zastosowano wspomagane komputerowo projektowanie leków, teorię funkcjonału gęstości, jakościową zależność struktura-aktywność i dynamikę molekularną w celu odkrycia silnych NNRTI.

Protokół

1. Szczegóły obliczeniowe — przygotowanie białka

  1. Kliknij ikonę Windows na ekranie monitora, a następnie wybierz Wszystkie aplikacje. Przewiń w dół do folderu Schrodinger, otwórz go i kliknij ikonę Maestro przedstawioną na rysunku 2B, a następnie wybierz otwórz, jak pokazano na rysunku 2B, aby uruchomić oprogramowanie.
  2. Pobierz wybraną strukturę białka, korzystając z karty Plik w oprogramowaniu. Z wyświetlonego krótkiego menu wybierz Pobierz PDB, jak pokazano na rysunku 3A, i wpisz wybrany kod PDB w polu tekstowym, jak pokazano na rysunku 3B. Kliknij przycisk pobierz, a wybrany plik PDB zostanie wyświetlony w oknie projektu.
  3. Alternatywnie można pobrać wybrane białko z Protein Data Bank na komputer lokalny, wpisując kod identyfikacyjny bazy białek (PDB ID) w polu wyszukiwania i klikając pobierz, aby zapisać plik PDB na komputerze lokalnym. Przejdź do karty Plik i wybierz opcję Importuj struktury. W interfejsie Importuj zlokalizuj pobrany plik PDB, jak pokazano na rysunku 3C, a następnie kliknij przycisk Importuj, zgodnie z oznaczeniem na rysunku 3D.
    UWAGA: Struktura białka zostanie otwarta jako struktura 3D w osobnym oknie, jak pokazano na rysunku 4.
  4. Przejdź do prawego górnego rogu oprogramowania, wybierz opcję  zadanie i wpisz przygotowanie białka w pasku wyszukiwania. Kliknij Protein Preparation Workflow w panelu po prawej stronie, jak pokazano na rysunku 5A.
  5. W wyświetlonym oknie Protein Preparation Workflow wpisz nazwę zadania jako nazwę pliku do zapisu i kliknij zielony przycisk Uruchom w prawym dolnym rogu, jak pokazano na rysunku 5B.
  6. Monitoruj postęp zadania, klikając przycisk Zadania w prawym górnym rogu, jak pokazano na rysunku 5C.
  7. Zaznacz przygotowane białko, kliknij na nie prawym przyciskiem myszy i wybierz split ligand, jak pokazano na rysunku 5D. Wybierz podział na ligandy, wodę i pozostałe elementy. Pozwoli to na uzyskanie komponentów białka jako niezależnych wpisów w nawigatorze obszaru roboczego 

2. Przygotowanie ligandu

  1. Pobrać pożądane związki chemiczne z bazy danych PubChem30, wpisując nazwę związku w polu wyszukiwania. Przejrzeć struktury i wybrać struktury trójwymiarowe (3D). Kliknąć Download w prawym górnym oknie, aby pobrać współrzędne struktury w formacie pliku danych strukturalnych (SDF). Jeśli struktura 3D jest niedostępna, pobrać strukturę 2D i użyć innych narzędzi do wygenerowania struktury 3D na podstawie struktury 2D.
  2. Kliknąć kartę  File w programie Schrodinger i wybrać Import Structures, jak pokazano na Rysunku 6A. Przejść do folderu, w którym zapisane są struktury w formacie SDF, aby załadować związki przeznaczone do przygotowania.
  3. Wybrać Task w prawym górnym rogu oprogramowania Schrodinger. Wpisać LigPrep w polu wyszukiwania i wybrać LigPrep w prawym oknie po lewej stronie, jak pokazano na Rysunku 6B.
  4. Wybrać Use Structures from , aby wskazać pliki z obszaru Workspace lub tabeli Project. Wybrać preferowane opcje w oknie LigPrep, zapisać plik na komputerze lokalnym i kliknąć run , aby uruchomić proces przygotowania ligandów, jak pokazano na Rysunku 6C.

3. Geometria i optymalizacja ligandów

  1. Uruchom oprogramowanie31 w celu optymalizacji geometrii pobranych struktur. Przejdź do karty File (Rysunek 7A) i wybierz Open, aby wybrać pobrany plik SDF z bazy danych PubChem.
    UWAGA: Plik zostanie załadowany w głównym oknie. Kliknij drugie fioletowe okno, aby zwizualizować te same związki chemiczne. Każde wprowadzone ustawienie wyświetli wynik w fioletowym oknie.
  2. Przejdź do karty Calculate tab i wybierz kartę Gaussian Calculation Setup Tab. Przejdź do karty Job type przedstawionej na Rysunku 7B i wybierz Optimization lub Opt+Freq. 
    UWAGA: W zależności od liczby (rozmiaru) atomów lub związków, najpierw przeprowadź optymalizację, a następnie obliczenia częstotliwości. W przypadku mniejszych związków wykonaj Opt+Freq. Im większa cząsteczka, tym więcej czasu zajmuje wykonanie Opt+Freq w porównaniu do optymalizacji, po której następuje obliczenie częstotliwości.
  3. Przejdź do karty Method i wybierz metody chemii kwantowej, a następnie z rozwijanych list w każdej sekcji wybierz odpowiedni globalny hybrydowy funkcjonał gęstości wymiany i korelacji Kohna-Shama , bazę funkcji, ładunek oraz spin.
  4. Przejdź do sekcji Title i określ nazwę badanego związku.
  5. Przejdź do karty Link 0 i określ preferowany limit pamięci (Memory limit) oraz liczbę procesorów współdzielonych (Shared processors). Odznacz pola Full path, jak pokazano na Rysunku 7C.
  6. Kliknij przycisk Edit na dole, aby zapisać plik wejściowy Gaussian, jak pokazano na Rysunku 3C. Zapisz plik w wybranej lokalizacji pod wybraną nazwą jako plik zadania Gaussian (gjf).
    UWAGA: Po zapisaniu pliku wejściowego Gaussian pojawi się okno pop-up wyświetlające zawartość pliku w Notatniku. Można edytować zawartość pliku, w tym zmienić ładunek i funkcjonał, jeśli opcje te były niedostępne w ustawieniach. Przygotowany plik zadania Gaussian zostanie użyty jako plik wejściowy do uruchomienia obliczeń optymalizacji i częstotliwości na komputerze lokalnym.
    Następnie przeprowadzono optymalizację geometryczną tych struktur z wykorzystaniem funkcjonału MN15-L32 i bazy funkcji 6-31++G(d,p)33. W przypadku wystąpienia problemów podczas procesu optymalizacji i obliczeń częstotliwości, zadanie można przesłać do wykonania z pomocą klastra obliczeniowego HPC.

4. Generowanie siatki receptora

  1. Przejdź do sekcji Tasks i wybierz receptor grid generation (generowanie siatki receptora)34. Interfejs generowania siatki receptora pokazany na Rycynie 8A służy do identyfikacji centrum aktywnego białka, w którym związany jest ligand krystaliczny. Kliknij Pick to identify the ligand (Wybierz, aby zidentyfikować ligand) i sprawdź, czy w strukturze obecny jest ligand kokrystalizowany, korzystając z wyskakującego powiadomienia u góry ekranu, jak pokazano na Rycynie 8B.
  2. Wybierz ligandy i/lub pozostałości (residues) w obszarze roboczym; wybrane związki będą wyróżnione kolorem niebieskim.
  3. Wybierz kartę ustawień w panelu generowania siatki receptora, aby zdefiniować grid box size (rozmiar pudełka siatki). Domyślny rozmiar pudełka siatki to 10 Å x 10 Å x 10 Å. Zmień wymiary pudełka w karcie Grid Box, wprowadzając żądane wartości bezpośrednio lub ręcznie dostosowując rozmiar pudełka w obszarze roboczym.
    UWAGA: W razie potrzeby można ustawić dodatkowe parametry w karcie  Advanced (Zaawansowane), jak pokazano na Rycynie 8C. Przejrzyj wszystkie ustawienia, aby upewnić się co do ich poprawności.
  4. Kliknij Run (Uruchom), aby rozpocząć generowanie siatki. Po zakończeniu procesu zapisz wygenerowane pliki siatki do późniejszych symulacji dokowania. Po zakończeniu zadania pojawi się komunikat powiadomienia, jak pokazano na Rycynie 8D.

5. Dokowanie molekularne

  1. Załaduj białko oraz przygotowane ligandy, przechodząc do Zadania, Dokowanie35or dokowanie ligandów (dokowanie Glide)), jak pokazano w Rycina 9A.
  2. Wczytaj plik siatki z kroku 4.1 powyżej i wybierz ligandy z obszaru roboczego, korzystając z opcji Use ligand from the opcja w Rysunek 9B.
  3. Wybierz preferowany metoda precyzji dokowania z ustawienia tabela przedstawiona w Rycina 9 C (domyślna precyzja dokowania to Standard Precision (SP)).
  4. Ustawić pole sił do OPLS4 w celu dokładnego modelowania oddziaływań molekularnych. Ustaw więzy (np. wiązania wodorowe) w Ograniczenia tab.
  5. Przejrzyj wszystkie ustawienia i zapisz zadanie lub plik dokowania. Kliknij URUCHOM** aby rozpocząć proces dokowania.
  6. Powiadomienie wskazuje, że zadanie zostało zakończone, a Tabela Projektu zostaje otwarta dla wyników dokowania. Przeanalizuj pozy pozycjonowania, wynikii oddziaływania.

6. Generowanie modelu 2D-QSAR

  1. Przejdź do strony centrum społeczności KNIME i wyszukaj AutoQsar w pasku wyszukiwania. Wybierz drugi wpis AutoQsar, który pojawia się na rysunku 10A.
  2. Kliknij ikonę pobierania przepływu pracy (workflow) i wybierz opcję pobierz przepływ pracy w menu wyskakującym widocznym na rysunku 10B. Zapisz przepływ pracy na komputerze lokalnym.
  3. Upewnij się, że program KNIME36 jest zainstalowany, a następnie uruchom go z komputera lokalnego. Przygotuj arkusz kalkulacyjny zawierający kanoniczne uśmiechy (SMILES), nazwę struktury, wartości aktywności/IC50 (w mikromolach) oraz -log(aktywności/IC) lub dowolny wybrany deskryptor chemiczny. Zapisz plik w formacie CSV (Comma-separated value) na komputerze lokalnym.
  4. Przejdź do menu Plik i zaimportuj przepływ pracy AutoQSAR37. W oknie wyskakującym kliknij Importuj przepływ pracy KNIME, aby wybrać pobrany przepływ pracy AutoQSAR, jak pokazano na rysunku 10C, a następnie wybierz Otwórz | Dalej | Zakończ.
  5. Przepływ pracy AutoQsar pojawi się w oknie Eksploratora KNIME w lewym górnym rogu. Kliknij dwukrotnie przepływ pracy, aby przenieść go do głównego okna.
  6. Kliknij dwukrotnie ikonę czytnika cząsteczek (do MAE), a następnie w karcie Ustawienia kliknij Dodaj plik(i), aby dodać arkusz z danymi wejściowymi lub parametrami.
  7. Wybierz kartę Polityka pamięci i kliknij Zapisz tabele na dysku | OK.
  8. Kliknij prawym przyciskiem myszy przycisk Wyodrębnij właściwości i wybierz Konfiguruj, jak pokazano na rysunku 10D.
  9. Wybierz i przefiltruj interesujące właściwości molekularne lub chemiczne, przenosząc je z okna wyklucz do okna dołącz. Kliknij Zastosuj.
  10. Kliknij prawym przyciskiem myszy element AutoQSAR Build Model i wybierz Konfiguruj. Wpisz nazwę modelu QSAR w oknie dialogowym.
  11. Kliknij przycisk ustawień i wybierz właściwość do dopasowania spośród tych dołączonych. Wybierz wartość losowego zbioru treningowego (np. 85%:15%) oraz liczbę modeli do zachowania. Kliknij Zastosuj.
  12. Kliknij prawym przyciskiem myszy dolny Czytnik cząsteczek i wybierz Konfiguruj. Załaduj ligandy testowe z przygotowanego pliku Excel CSV zawierającego ligandy do przetestowania.
  13. Kliknij prawym przyciskiem myszy Wyodrębnij właściwości i ustaw parametry według własnego uznania.
  14. Kliknij prawym przyciskiem myszy górny Czytnik cząsteczek i wybierz wykonaj. 
    UWAGA: Pomarańczowa kropka oznacza rozpoczęcie procesu, a zielone światło wskazuje na pomyślne wykonanie procesu. Kolejne kroki będą realizowane, aż wszystkie węzły zostaną pomyślnie wykonane.
  15. Wykonaj drugi Czytnik cząsteczek dla ligandów testowych.
  16. Zapisz folder zip z wynikami ze zbioru testowego i treningowego po  pomyślnym wykonaniu wszystkich węzłów.
    UWAGA: Wybierz najlepiej działający model na podstawie wyników walidacji krzyżowej, takich jak SD – odchylenie standardowe, R2 – korelacja zbioru treningowego między rzeczywistymi a przewidywanymi wartościami aktywności oraz Q2 – korelacja rzeczywistej i przewidywanej aktywności zbioru testowego.
  17. Oceń wydajność przepływu pracy, analizując sprawność modelu za pomocą  węzłów oceniających (scorer) lub statystycznych.
  18. Aby zinterpretować wyniki, zidentyfikuj główne lub kluczowe deskryptory, korzystając z istotności cech (feature importance). Zwizualizuj wyniki w formie wykresu rozrzutu lub wykresu słupkowego, aby określić wydajność i relację deskryptorów. Zapisz przepływ pracy i wyeksportuj dane wynikowe, przewidywania oraz wizualizacje w formie wykresów rozrzutu i/lub plików CSV.

7. Wyliczanie

  1. Przejdź do sekcji Task i wyszukaj program Ligand Designer, jak pokazano na Rysunku 11A.
  2. W nawigatorze obszaru roboczego wybierz parę złożoną z zadokowanej białka i liganda. W oknie Ligand Designer kliknij Analyze Workspace. Aby wygenerować i ocenić nowe ligandy, z listy pojawiających się przepływów pracy (workflows) wybierz Isostere scanning, jak pokazano na Rysunku 11B , co oznacza zastosowanie metody wzrostu (growing method), która rozbudowuje ligand poprzez dodawanie fragmentów do istniejących części cząsteczki.
  3. Po ustawieniu opcji enumeracji, kliknij enumerate w pojawiającym się oknie powiadomień Isostere Scanning. Powtórz krok 6.2 dla tego samego białka i innego liganda, aż wszystkie ligandy zostaną poddane temu samemu procesowi. Przejrzyj wygenerowane ligandy w Tabeli Projektów (Projects Table).
  4. Zapisz struktury, eksportując zaprojektowane ligandy.
    UWAGA: Metoda enumeracji wygenerowała zestaw wyników dokowania (docking scores) po zakończeniu symulacji.
  5. Indywidualnie wybierz zenumerowane związki o bardziej ujemnych wynikach dokowania i przeprowadź ponowne dokowanie, stosując procedurę opisaną w sekcji 4.

8. HOMO i LUMO

  1. Otwórz oprogramowanie i prześlij zoptymalizowany ligand w formacie pliku zadania Gaussian.
  2. Przejdź do narzędzia i wybierz Edytor MO w czerwonych i zielonych kropkachjak pokazano w Rysunek 12A.
    UWAGA: Metoda enumeracji wygenerowała zestaw wyników dokowania po zakończeniu symulacji.
  3. W edytor MO panel poniżej metoda, załaduj plik FChk poprzez kliknięcie w załaduj Mos z istniejącego pliku Chk lub FChkjak pokazano w Rycina 12B.
  4. Kliknij w wizualizować zakładka aktualizacjaOdczekać ok. 10 s do momentu pojawienia się powierzchni prądu, jak pokazano w Rycina 12D.
  5. Kliknij jeden z dwóch pola wyboru obok podświetlonych na żółto cyfr, aby wybrać jedną z nich HOMO lub LUMO38 powierzchnia do wyświetlenia w sąsiednim oknie.
  6. Kliknij prawym przyciskiem myszy na fioletowym tle i wybierz PlikKliknij w zapisz plik obrazu aby zapisać obraz aktualnej powierzchni, jak pokazano w Rycina 12E.
  7. Kliknij prawym przyciskiem myszy na fioletowym tle i wybierz Wyświetl Wybierz format wyświetlania aby zmienić tło w Ogólne zakładka, przezroczystość powierzchni w Powierzchnia zakładka, rozmiar czcionki i kolor w karcie Text, jakość obrazu oraz preferowany układ powierzchni w karcie Molecule, jak pokazano w Rysunek 12F.
  8. Alternatywnie można wygenerować plik cube i wgrać plik FChk w programie GaussView. Przejdź do Wyniki karta i wybór Powierzchnie/KonturyKliknij w działania z kostką i załaduj sześcian. Kliknij w Działania powierzchniowe i wybierz nową powierzchnię. Powtórz krok 8.7, aby edytować powierzchnię.
    UWAGA: Im mniejsza przerwa energetyczna między HOMO a LUMO39 (różnica między LUMO a HOMO), tym bardziej reaktywną powinna być dana cząsteczka. Oblicz przerwę energetyczną (Eprzerwa) przy użyciu równania 1 dla każdej cząsteczki.
    figure-protocol-1

9. Symulacja dynamiki molekularnej — przygotowanie i minimalizacja układu

  1. Upewnij się, że pakiet Schrodinger Suite został zainstalowany na lokalnym komputerze, a następnie załaduj struktury kompleksów białko-ligand do obszaru roboczego.40.
  2. Kliknij w Zadanie przycisk i wybierz Konfigurator Systemu Desmond41W w ramach konfigurator systemu panel, wybierz Solwatacja przejdź do karty i wybierz Predefiniowany model rozpuszczalnika jak pokazano w Rycina 13A, Kształt pudełka jak pokazano w Rysunek 13Boraz rozmiar pudełka metoda obliczeń przedstawiona w Rycina 13C42odpowiadające kompleksowi białko-ligand.
  3. Wybierz jony przejdź do karty i kliknij przelicz w celu zneutralizowania układu poprzez dodanie przeciwjonów i ustawienie pożądane stężenie soli.
  4. Wybierz OPLS43,44 wybranym jako pole sił.
  5. Nadaj zadaniu odpowiednią nazwę i zapisz plik zadania na komputerze lokalnym. Kliknij uruchom aby przekazać zadanie do przygotowania.
    UWAGA: Należy upewnić się, że nazwa zadania została zapisana. Należy użyć szczegółowej nazwy.
  6. Równoważenie i produkcja
  7. Wyświetl projekt w obszarze roboczym po przygotowaniu systemu. Wybierz kompleks białko-ligand z Nawigatora obszaru roboczego (Workspace Navigator), przejdź do zadanie, i wybierz dynamika molekularna (Desmond).
  8. Wczytaj kompleks ligand-białko z obszaru roboczego w panelu dynamiki molekularnej. Wybierz żądane oś czasu symulacji z symulacja zakładka Wybierz NPT jako klasa zespołowa45.
  9. Nazwij zadanie w odpowiedni sposób i opisz je. Kliknij w Zamknij aby wyjść z okna dynamiki molekularnej.
  10. Prześlij przygotowane zlecenie symulacji dynamiki molekularnej za pomocą lokalnego terminala. Po zakończeniu otwórz ukończone zlecenie i kontynuuj symulację, wydłużając czas z początkowo ustawionego przedziału do pożądanego czasu symulacji.45na przykład 100 ns, 200 ns.
    UWAGA: Powtórz krok 9 dla pozostałych kompleksów białko-ligand.
    1. Otwórz zakończone zadanie w oprogramowaniu Schrodinger Maestro i połącz różne ramy czasowe symulacji, jeśli uruchomiono oddzielne zadania. Przejdź do ZADANIE przycisk i wybierz schemat interakcji symulacyjnychWczytaj połączone pliki lub pojedynczy plik, aby zwizualizować trajektorię.
  11. Wygeneruj raport zawierający wizualizacje oraz inne szczegółowe wyniki wyjściowe.

10. Mechanika molekularna z uogólnionym modelem Borna i powierzchnią dostępną dla rozpuszczalnika (MM-GBSA)

  1. Otwórz plik trajektorii i odtwórz trajektorię. Zwizualizuj, gdzie kompleks białko-ligand uległ equilibracji, i zanotuj liczbę klatek. Prześlij zadanie do przetworzenia za pomocą terminala.
  2. Powtórz sekcję 9 (przygotowanie symulacji dynamiki molekularnej białka) dla drugiego kompleksu białko-ligand.
  3. Połącz lub przejrzyj zawartość pliku wyjściowego, aby przeanalizować wygenerowane wyniki. Odczytaj wartość ΔG Average, aby uzyskać energię swobodną wiązania kompleksu białko-ligand.
  4. Pobierz plik CSV, aby zwizualizować wkłady różnych cząsteczek wewnątrzmolekularnych.
  5. Aby obliczyć swobodną energię wiązania kompleksu, uwzględnij poszczególne parametry termodynamiczne i desolwatacyjne, w tym energię wiązania (ΔGbind), kulombowski model solwatacji (ΔGbind Coulumb), niepolarny składnik solwatacji (ΔGbind Lipo), korektę wiązań wodorowych (ΔGbind Hbond), wiązanie kowalencyjne (ΔGbind Covalent), korektę upakowania π-π (ΔGbind Packing), energię solwatacji elektrostatycznej Generalized Born (ΔGbind sol GB) oraz oddziaływanie van der Waalsa (ΔGbind vdW).
    1. Wyznacz końcową wartość ΔGbind poprzez uśrednienie wartości ΔGbind określonych dla każdego snapshotu w symulacji MD, zgodnie z równaniem 2.

figure-protocol-2

Wyniki

Generowanie siatki receptora i dokowanie molekularne

Do właściwej charakterystyki miejsca wiązania przed kolejnym etapem dokowania wykorzystano narzędzie do generowania siatki receptorów w programie Maestro. Do zdefiniowania siatki użyto ligandu kokrystalizowanego. W dokowaniu molekularnym zastosowano ustawienie Glide o precyzji SP. Do przygotowania ligandów do dokowania z wykorzystaniem pola siłowego OPLS4 użyto narzędzia LigPrep w programie Schrödinger Maestro. W symulacjach dynamiki molekularnej w programie Desmond zastosowano pole siłowe OPLS4. Pola siłowe są fundamentalne dla klasycznych symulacji molekularnych, a ich dokładność ma kluczowe znaczenie dla jakości symulacji wiązania białko-ligand w procesie odkrywania leków. W przypadku OPLS4 przypisanie ładunków i parametrów wykonano w programie Schrödinger Maestro. Zastosowanie parametrów OPLS4 doprowadziło do znaczącej poprawy zarówno w porównaniach energetycznych, jak i geometrycznych w odniesieniu do parametrów domyślnych OPLS2005.

Wiązało się to z dokowaniem specyficznych ligandów znanych z powinowactwa do białka docelowego HIV-1, co pozwoliło na dokładną analizę oddziaływań między cząsteczkami ligandów a resztami receptora. Tabela 1 przedstawia podsumowanie klasyfikacji związków dla modelowania QSAR.

Po dokonaniu dokowania HBY561, protokół dokowania oceniono poprzez porównanie redokowanego ligandu z ligandem znajdującym się w centrum aktywnym krystalicznego białka 1HQU. Struktury zadokowanego i skrystalizowanego HBY561 znajdują się w Pliku uzupełniającym 1 (Rycina uzupełniająca S1 oraz Rycina uzupełniająca S2). Aby ocenić podobieństwo między pozami dokowania a strukturami referencyjnymi, powszechnie przyjmuje się, że wartość średniego kwadratu odchylenia (RMSD) mniejsza niż 2.0 Å stanowi kryterium wiarygodnych wyników dokowania. Próg ten wskazuje, że przewidywana struktura ściśle pokrywa się z danymi eksperymentalnymi. W niniejszym badaniu ligandy wykazały wartość RMSD wynoszącą 1.27 Å przy porównaniu struktury referencyjnej z pozą dokowania, co zilustrowano kolorem czerwonym na Rycinie uzupełniającej S3. Wykazało to, że protokół dokowania był wystarczający dla niniejszej pracy, w związku z czym wszystkie ligandy zadokowano przy użyciu tych samych ustawień. Analizując wyniki dokowania przedstawione w Tabeli 2, Efavirenz i Etravirina wykazały najkorzystniejsze wyniki wynoszące odpowiednio -10.432 eV i -9.647 eV w stosunku do wyniku dokowania ligandu współkrystalizowanego HBY561 (-9.242 eV). Rycina uzupełniająca S4 z Pliku uzupełniającego 1 przedstawia diagramy oddziaływań ligandów między oryginalnym białkiem HIV-1 a ligandem krystalicznym, Etraviriną oraz Efavirenzem.

Wiązania wodorowe, stosowanie π-π oraz oddziaływania hydrofobowe są głównymi siłami przyczyniającymi się do wiązania. Konkretny przykład wiązania wodorowego pojawił się pomiędzy HBY561 a wyznaczonym białkiem 1HQU, z wyraźnym udziałem reszty aminokwasowej LYS101. Ten wzór wiązań odzwierciedlał obserwacje poczynione w przypadku ligandów Efavirenz i Etravirine, jak pokazano na Rysunku 14.

Ponadto oddziaływania hydrofobowe były kluczowe dla wiązania w kilku miejscach białka z udziałem HBY561, etrawiryny i efawirenzu, obok sił międzycząsteczkowych, wiązań wodorowych i stosowania π-π. Stosowanie π-π zaobserwowano pomiędzy TYR318 a pierścieniem aromatycznym w efawirenzu. Zarówno wiązania wodorowe, jak i stosowanie π-π były niezbędne do utrzymania połączeń wiążących pomiędzy ligandami a białkiem. Diagramy oddziaływań ligandów HBY_561, newirapiny, dorawiryny, efawirenzu i etrawiryny przedstawiono w Supplemental File 1-SupplementalFigureS5.

Te oddziaływania międzycząsteczkowe wpływają na interakcje białko-ligand i mają kluczowe znaczenie dla opracowywania leków zmniejszających AMR w przypadku HIV-1. Ich rola w zwiększaniu powinowactwa wiązania, swoistości oraz mechanizmu działania pomaga w projektowaniu leków, które mogą skutecznie celować w wirusa i go hamować, rozwiązując tym samym rosnący problem oporności na środki przeciwdrobnoustrojowe w kontekście leczenia HIV-1.

Przygotowanie zbioru danych 2D-QSAR

Fazy uczenia i testowania objęły 94 związki. Związki te podzielono na cztery klasy, z których każda reprezentowała ligandy związane ze specyficznym badanym białkiem. W procesie uczenia wykorzystano przepływ pracy KNIME AutoQSAR, w którym jako trzy deskryptory do niniejszego badania wybrano aktywność, HOMO oraz LUMO.

Generowanie 2D-QSAR

Wartości aktywności wszystkich 94 cząsteczek zostały określone na podstawie danych eksperymentalnych (Tabela uzupełniająca S1). Wyniki uzyskane z modelu QSAR znajdują się w Tabeli 3. Do modelowania QSAR wykorzystano związki z Klasy 1 wymienione w Tabeli uzupełniającej S2, z uwzględnieniem dołączonych grup Ar oraz ich odpowiednich wartości aktywności. Wyniki dla czterech klas wskazują, że najwyższy wynik 0,8223, R2 na poziomie 0,815 oraz Q 2 wynoszące 0,8182 osiągnięto dla Klasy 1. Jest to zgodne z wcześniejszymi kryteriami dążenia do uzyskania R2 bliskiego 1 oraz Q2 większego niż 0,746. W związku z tym do trenowania naszego modelu QSAR wybrano klasę 1. Chociaż modele dla klas 3 i 4 wykazały doskonałe wartości korelacji R2 wynoszące odpowiednio 0,8172 i 0,6673, nie dorównywały one wydajności Klasy 1.

W celu dalszej walidacji stabilności proponowanego modelu QSAR przeprowadzono korelację krzyżową zgodnie z kryterium, według którego różnica między wynikiem Q2 (0,8223) a R2 powinna być mniejsza lub równa 0,347. W przypadku naszego proponowanego modelu dla klasy 1 różnica ta wynosi 0,0038. Wykres rozrzutu przedstawiający aktywność obserwowaną w stosunku do aktywności przewidywanej znajduje się na Rysunku 15.

Przerwa energetyczna HOMO-LUMO

Wyznaczenie przerwy energetycznej między najniższym nieobsadzonym orbitalem molekularnym a najwyższym obsadzonym orbitalem molekularnym, powszechnie znanej jako przerwa energetyczna HOMO-LUMO, odgrywa kluczową rolę w charakteryzacji reaktywności chemicznej i stabilności kinetycznej cząsteczki w kontekście sześciu związków NNRTIs. Orbitale molekularne graniczne odgrywają kluczową rolę w ułatwianiu oddziaływań przeniesienia ładunku z miejscem wiązania białka HIV. Potwierdzenie minimum energii zapewniono poprzez analizę częstotliwości wibracyjnych i potwierdzenie braku częstotliwości ujemnych lub urojonych; następnie dla każdego minimum energii wyznaczono wartości HOMO i LUMO. Wyższa wartość HOMO oznacza większą zdolność cząsteczki do bycia donorem elektronów, podczas gdy niższa wartość sugeruje, że działa ona jako słaby akceptor elektronów. Rycina dodatkowa S6 przedstawia przerwy energetyczne HOMO_LUMO dla sześciu zoptymalizowanych związków NNRTIs. Ponadto zmniejszona przerwa energetyczna między poziomami energetycznymi HOMO i LUMO silnie wpływa na międzycząsteczkowe oddziaływania przeniesienia ładunku, które występują między badanymi cząsteczkami ze względu na ich silną zdolność do akceptowania elektronów oraz bioaktywność48.

Trend wartości przerwy energetycznej, przedstawiony w tabeli 4, wykazuje kolejność malejącą: Efavirenz > Etravirine > HBY-561 > Nevirapine > Delavirdine > Doravirine > Rilpivirine. Znaczna przerwa energetyczna zaobserwowana dla Efavirenz i Etravirine oznacza, że analizy wyników dokowania ujawniają korelację między bioaktywnością a przerwą HOMO-LUMO. Co istotne, potencjał przeciwwirusowy wzrasta wraz z większymi wartościami przerwy HOMO-LUMO. Wskazuje to nie tylko na stabilność związków, ale także na ich zdolność do tworzenia trwałych oddziaływań z receptorem. Przerwa HOMO-LUMO odgrywa istotną rolę w zrozumieniu bioaktywności cząsteczek, szczególnie w kontekście projektowania leków przeciwko HIV-1.

Liczenie

Do wyliczania bibliotek wirtualnych stworzono różne narzędzia. Wśród narzędzi wykorzystywanych do wyliczania znajduje się oprogramowanie Schrödinger. Opiera się ono na metodzie przeskakiwania rdzenia (core hopping), w której biblioteki są tworzone poprzez zastępowanie jednej lub kilku grup przyłączonych do struktury rdzeniowej fragmentami związków odczynnikowych49. Narzędzie do wyliczania w programie Maestro v13.1 zostało użyte do dodania niestandardowych grup bocznych lub atomów do każdego z sześciu NNRTIs. Nowe związki zostały również wykorzystane do przewidywania aktywności. W wartościach przewidywanej aktywności wyliczonych cząsteczek odnotowano poprawę w porównaniu z początkowo zoptymalizowanymi cząsteczkami NNRTIs, co przedstawiono w Table 5.

Poprawa wartości aktywności wykazana dla ligandów z enumeracji była wynikiem procesu enumeracji, ponieważ dodana niestandardowa grupa R wpłynęła na siły oddziaływania nowo zaproponowanego związku z oryginalnym białkiem. Związki z enumeracji przygotowano do obliczeń mechaniki kwantowej poprzez optymalizację tych cząsteczek i obliczenie ich częstotliwości wibracyjnych. Obliczono ich przerwy energetyczne i porównano je z przerwami energetycznymi NNRTIs. Ogólna obserwacja, przedstawiona w Tabeli 6, wskazuje, że związki z enumeracji są bardziej stabilne niż ich zoptymalizowane odpowiedniki.

W Tabeli 6 zaobserwowano, że przerwa energetyczna wyliczonych związków w porównaniu ze związkami zoptymalizowanymi wykazała podobną tendencję. Wskazuje to na to, że właściwości chemiczne wyliczonych cząsteczek pozostały niezmienione, niezależnie od rotacji konformacyjnej. Dzięki temu były one w stanie utrzymać istotne oddziaływania międzycząsteczkowe z resztami aminokwasowymi białka.

Przed przeprowadzeniem procesu enumeracji dla NNRTIs, zgodnie z opisem w sekcji 6 protokołu, zaobserwowane wyniki wyjściowe posiadały początkowe wyniki dokowania z procesu enumeracji, przedstawione w Tabeli 7 w kolumnie z wynikami dokowania enumerowanych związków. Jak wyjaśniono w sekcji 4 protokołu, proces dokowania molekularnego został przeprowadzony w celu walidacji proponowanego wyniku dokowania dla enumerowanych związków. Zaobserwowano, że po ponownym dokowaniu enumerowanych związków, nowe wyniki dokowania uległy poprawie, co pokazano w kolumnie „redocked enumerated ligands”. Wyniki ponownego dokowania dla enumerowanych związków porównano z wynikami dokowania oryginalnych NNRTIs przedstawionymi w kolumnie „original docking score” w Tabeli 7. Zaobserwowano, że w przypadku związków enumerowanych, wyniki dokowania dla HBY_561, Etrawiryny, Efavirenzu i Dorawiryny były lepsze niż dla ich odpowiednich związków zoptymalizowanych. Jednak Delawirydyna posiada ten sam wynik dokowania zarówno w wersji enumerowanej, jak i zoptymalizowanej.

Dynamika molekularna

HBY 561 tworzy trzy wiązania wodorowe z LYS101, silnie elektroujemnymi atomami N oraz grupą OH. Ligand krystaliczny, czyli trzecia cząsteczka, tworzy jednocześnie dwa wiązania wodorowe z siarką i wodorem oraz trzecie wiązanie wodorowe z GLU138. Co istotne, istotny wkład w siły międzycząsteczkowe mają oddziaływania π-π oraz oddziaływania hydrofobowe. Ponadto obecność dodatkowych wiązań wodorowych ma kluczowe znaczenie dla dynamiki molekularnej, co znajduje odzwierciedlenie w wynikach na wykresach średniokwadratowego odchylenia (RMSD). Łącznie zaobserwowano, że efawirenz tworzy trzy wiązania wodorowe: z bardzo elektroujemnym atomem azotu, atomem tlenu w drugim pierścieniu niearomatycznym oraz z pierścieniem benzenowym i TYR318 pomiędzy silnie elektroujemnym N w pierścieniu centralnym liganda. Widoczne jest drugie wiązanie wodorowe pomiędzy N z pierścienia a grupą OH i LY101. Etrawiryna wykazuje trzy wiązania wodorowe z LYS101. Dorawiryna tworzy wiązanie wodorowe z GLU138. Newirapina wykazuje dwa wiązania wodorowe z LY101.

W tym przypadku pozostałością aminokwasową, z którą oddziałują wszystkie ligandy, jest LYS101. Mimo że ich struktury się różnią, wszystkie one oddziałują z tą samą pozostałością aminokwasową. Symulacje MD przeprowadzono zgodnie z parametrami wymienionymi w sekcji 2.8 protokołu, aby określić, jak silnie lub słabo każdy ligand (NNRTI oraz wymienione NNRTI) wiąże się z miejscem aktywnym 1HQU. Oddziaływania ligandów przedstawione na Rysunku 16 wskazują na silne wiązania wodorowe pomiędzy aminokwasem LY101 białka a cząsteczkami HBY 561, Newirapiną, Efavirenzem i Etrawiryną. Jak wynika z porównania w Tabeli 7, te silne oddziaływania odpowiadają za wysokie wyniki dokowania każdego z związków.

W celu oceny skuteczności wiązania każdego liganda, w tym NNRTI oraz wyliczonych NNRTI, w centrum aktywnym 1HQU, przeprowadzono symulacje dynamiki molekularnej (MDS). W szczególności wybrano cztery wyliczone związki, które wykazały lepsze wyniki dokowania niż oryginalne zoptymalizowane kombinacje. Wybrane związki poddano MDS jako metodzie walidacji, aby zbadać i zaobserwować reakcję każdej cząsteczki z białkiem HIV-1 w określonym czasie, uwzględniając oddziaływania międzyatomowe w obecności liganda.

Przed rozpoczęciem MDS dla niedawno wyliczonych glikanów konieczne było potwierdzenie przydatności protokołu symulacji dla naszego układu. W tym celu pierwszym krokiem było przeprowadzenie MDS wolnego białka 1HQU. W centrum aktywnym białka nie znajdował się żaden ligand (Rysunek 17A oraz Rysunek dodatkowy S7). Do około 60 ns występują znaczne fluktuacje w strukturze białka, powodujące przesunięcia RMSD dla Cα wynoszące do 4,5 Å; po tym czasie białko wydaje się stabilizować, z fluktuacjami RMSD na poziomie ~3,5 Å do 200 ns. Stabilizacja ta przekonała nas, że protokół MD będzie odpowiedni dla naszych kompleksów białko-ligand, które zostaną przeanalizowane w następnej sekcji.

Symulacje MDS trwające 200 ns dla Etraviryny i wyliczonej Etraviryny (Rycina 17B,C) wykazały fluktuacje średniokwadratowego odchylenia (RMSD) Etraviryny bliskie 5.0 Å oraz stabilizację na poziomie 4.5 Å. Wyliczona Etraviryna wykazała fluktuacje RMSD wynoszące 4.5 Å i stabilizację na poziomie 3.5 Å. Stabilizacja ta wskazuje, że wyliczona Etraviryna może być potencjalnym ligandem NNRTI w leczeniu HIV/AIDS. Trajektoria z RMSD poniżej 5 Å oznacza silny efekt wiązania między białkiem miejsca aktywnego a ligandem. Obserwacja ta była właściwa dla wszystkich wcześniej wspomnianych związków, z wyjątkiem Newirapiny i Dorawiryny spośród związków wyliczonych (Plik uzupełniający 1: Rycina uzupełniająca S8, Rycina uzupełniająca S9, Rycina uzupełniająca S10 oraz Rycina uzupełniająca S11).

Rysunek 18 oraz Rysunek 19 stanowią dalszą analizę wiązania etrawiryny, wyliczonej etrawiryny oraz białka. Analizowane dane obejmują histogram interakcji oraz kontaktów ligand-białko i białko-ligand. Histogram kontaktów interakcyjnych dla każdego z ligandów związanych z białkiem jest bezpośrednio powiązany z odpowiadającymi im siłami oddziaływania między resztami aminokwasowymi białka a ligandem. Wysoka liczność LYS101 jest bardzo wyraźna w przypadku etrawiryny i wyliczonej etrawiryny, co objawiło się widocznym, grubym pomarańczowym pasem. Słabo dostrzegalny jasnopomarańczowy pas, znajdujący się w dolnej części wykresu dla wyliczonej etrawiryny, zaobserwowano w korelacji z TYR181. Zależność ta wskazuje na istnienie dwóch międzycząsteczkowych sił przyciągania w obrębie GLU138. Te pozytywne obserwacje dla ligandów oraz ich form wyliczonych zostały porównane w celu wyłonienia lepszego liganda jako potencjalnego związku NNRTI. Na podstawie przedstawionych wyników, wyliczona etrawiryna wykazuje potencjał do zastosowania w leczeniu HIV/AIDS. 

Obliczenia mechaniki molekularnej z uogólnionym modelem Borna i powierzchnią dostępną dla rozpuszczalnika (MM-GBSA)

W niniejszym badaniu głównym źródłem wkładu energetycznego do wolnej energii wiązania, ΔGbind, był udział oddziaływań van der Waalsa, ΔGVdW, . Wymieniona etrawiryna wykazuje wyższą wartość ΔGVdW wynoszącą -66,146 kcal/mol w porównaniu do znanego odpowiednika NNRTI, którego ΔGVdW wynosi -64,669 kcal/mol. Wyższa wartość ΔGHbond wynosząca -2,541 kcal/mol dla wymienionej etrawiryny wskazuje na istotny udział sił przyciągania wodorowego pomiędzy ligandem a białkiem. Wkłady ΔGCoulomb oraz ΔGCovalent dla wymienionej etrawiryny (-17,976 oraz 2,807 kcal/mol) były znacznie większe niż w przypadku znanego odpowiednika NNRTI, dla którego wartości te wynosiły odpowiednio -11,196 i 2,491.

Wyniki obserwowane podczas wiązania Etravirine oraz wyliczonego Etravirine z białkiem 1HQU są w pewnym stopniu spójne. Stwierdzono, że wyliczony Etravirine jest lepszy od swojego odpowiednika, Etravirine, ze względu na dwa dodatkowe wiązania wodorowe. Badanie wykazało również, że wyliczony Etravirine jest preferowany w wiązaniu z zidentyfikowaną kieszenią. Bardziej ujemna wartość ΔGbind (-89,684 kcal/mol) dla wyliczonego Etravirine w porównaniu do Etravirine (-80,551 kcal/mol) wskazuje, że wyliczony Etravirine jest dobrym inhibitorem RT HIV-1.

figure-results-1
Rycina 1: Struktury chemiczne sześciu zatwierdzonych przez amerykańską Agencję ds. Żywności i Leków (FDA) nienukleozydowych inhibitorów odwrotnej transkryptazy stosowanych w leczeniu wirusa niedoboru odporności typu 1 (HIV-1). NVP = Newirapina; DLV = Delawirdyna; EFV = Efavirenz; ETV = Etrawiryna; RPV = Rilpiwiryna; DOR = Dorawiryna. Prosimy kliknąć tutaj, aby wyświetlić powiększoną wersję tej ryciny.

figure-results-2
Rycina 2: Uruchamianie aplikacji Maestro Schrodinger w systemie Windows na komputerze lokalnym. (A) Lokalizacja aplikacji Maestro Schrodinger na komputerze lokalnym. (B) Sposób otwierania i uruchamiania aplikacji Maestro Schrodinger na komputerze lokalnym. Kliknij tutaj, aby wyświetlić powiększoną wersję tej ryciny.

figure-results-3
Rycina 3: Importowanie struktury z pliku PDB z lokalnego komputera do okna projektu w programie Schrödinger. (A) Funkcja importowania struktury w Maestro Schrodinger. (B) Pole tekstowe PDB ID. (C) Pobrany plik PDB na lokalnym komputerze. (D) Przycisk importowania umożliwiający zaimportowanie wprowadzonego pliku PDB ID. Kliknij tutaj, aby wyświetlić powiększoną wersję tej ryciny.

figure-results-4
Rysunek 4: Struktura pliku PDB zaimportowanego do okna projektu Schrodinger. Kliknij tutaj, aby wyświetlić powiększoną wersję tego rysunku.

figure-results-5
Rysunek 5: Schemat przygotowania białka. (A) Interfejs wyszukiwania w module Protein Preparation Workflow. (B) Zapisywanie nazwy pliku zadania i inicjowanie procesu przygotowania białka. (C) Okno monitorowania uruchomionych zadań. (D) Rozbicie ligandu na poszczególne komponenty. Kliknij tutaj, aby wyświetlić powiększoną wersję tego rysunku.

figure-results-6
Rycina 6: Schemat przygotowania ligandów. (A) Importowanie struktur z lokalnego komputera do okna projektu Protein preparation w programie Schrodinger. (B) Wyszukiwanie procesu Ligand preparation. (C) Okno schematu przygotowania ligandów. Prosimy kliknąć tutaj, aby wyświetlić powiększoną wersję tej ryciny.

figure-results-7
Rysunek 7: Geometria i schemat optymalizacji. (A) Okno menu GaussView do optymalizacji geometrii. (B) Dostępne Typy Zadań (Job Types) w zakładce Oblicz (Calculate) programu GaussView. (C) Opcje dostępne w zakładce Link0 w programie GaussView. Kliknij tutaj, aby wyświetlić większą wersję tego rysunku.

figure-results-8
Rysunek 8: Schemat generowania siatki Glide oraz workflow dokowania molekularnego. (A) Interfejs workflow generowania siatki receptora. (B) Powiadomienie wyskakujące dotyczące wyboru atomu wewnątrz ligandu. (C) Zaawansowane ustawienia receptora. (D) Powiadomienie o zakończeniu zadania. Kliknij tutaj, aby wyświetlić powiększoną wersję tego rysunku.

figure-results-9
Rysunek 9: Dokowanie ligandów w programie Glide. (A) Interfejs wyszukiwania dokowania ligandów Glide. (B) Interfejs dokowania ligandów. (C) Ustawienia precyzji dla dokowania Glide. Kliknij tutaj, aby wyświetlić powiększoną wersję tego rysunku.

figure-results-10
Rycina 10: Przepływ pracy przygotowania QSAR w programie KNIME. (A) Wyszukiwanie węzła AutoQSAR na stronie KNIME community hub. (B) Przycisk pobierania w celu uzyskania węzła AutoQSAR KNIME. (C) Importowanie pobranego przepływu pracy AutoQSAR KNIME. (D) Ustawienia konfiguracji dla ligandów podczas budowania modelu QSAR. Kliknij tutaj, aby wyświetlić większą wersję tej ryciny.

figure-results-11
Rysunek 11: Wyliczanie ligandów przy użyciu modułu Ligand Designer w programie Maestro Schrodinger. (A) Wyszukiwanie opcji Ligand Designer w programie Schrodinger. (B) Lista przepływów pracy (workflows) dla procesu wyliczania. Kliknij tutaj, aby zobaczyć powiększoną wersję tego rysunku.

figure-results-12
Rycina 12: Schemat generowania HOMO-LUMO. (A) Dostęp do opcji edytora cząsteczkowego w karcie Tools w programie GaussView. (B) Wczytywanie istniejącego pliku Chk lub FChk w celu wygenerowania granicznych orbitali molekularnych. (C) Ilustracja granicznych orbitali HOMO i LUMO w programie GaussView. (D) Okno wizualizacji służące do wyświetlania granicznych orbitali HOMO i LUMO. (E) Zapisywanie granicznych orbitali HOMO i LUMO. (F) Interfejs formatu wyświetlania w programie GaussView. Kliknij tutaj, aby wyświetlić powiększoną wersję tej ryciny.

figure-results-13
Rycina 13: Schemat przygotowania, konfiguracji i wykonania symulacji dynamiki molekularnej. (A) Desmond System Builder: Opcje solwatacji służące do wyboru sztywnych modeli wody. (B) Desmond System Builder: Opcje granic (Boundary) służące do wyboru kształtu pudełka. (C) Desmond System Builder: Metoda obliczania rozmiaru pudełka. Kliknij tutaj, aby zobaczyć powiększoną wersję tej ryciny.

figure-results-14
Rycina 14: Schematy oddziaływań ligandów między 1HQU a 3 najlepiej dopasowanymi ligandami. Siły oddziaływań między białkiem (1HQU) a (A) ligandem krystalicznym (HBY561), (B) efawirenzem oraz (C) etrawiryną. Siły te wpływają na oddziaływania białko-ligand i mają kluczowe znaczenie dla opracowywania leków zmniejszających oporność przeciwdrobnoustrojową w przypadku HIV-1. Ich rola w zwiększaniu powinowactwa wiązania, specyficzności i mechanizmu działania pomaga w projektowaniu leków, które mogą skutecznie celować w wirusa i go hamować, rozwiązując tym samym narastający problem oporności przeciwdrobnoustrojowej w kontekście leczenia HIV-1. Kliknij tutaj, aby wyświetlić powiększoną wersję tej ryciny.

figure-results-15
Rycina 15: Wykres punktowy przedstawiający aktywność obserwowaną w stosunku do aktywności przewidywanej dla klasy 1 modelu QSAR. Wykres przedstawia dopasowanie klasy 1 jako zbioru treningowego i związków NNRTI jako zbioru testowego w celu uzyskania przewidywanej wartości aktywności. Skróty: NNRTI = nie nukleozydowe inhibitory odwrotnej transkryptazy; QSAR = ilościowa zależność struktura-aktywność. Kliknij tutaj, aby wyświetlić większą wersję tej ryciny.

figure-results-16
Rycina 16: Schematy oddziaływań ligandów. Siły oddziaływań między białkiem a (A) ponumerowanym ligandem krystalicznym HBY_561, (B) ponumerowaną Newirapiną, (C) ponumerowaną Dorawiryną, (D) ponumerowanym Efavirenzem oraz (E) ponumerowaną Etrawiryną. Kliknij tutaj, aby wyświetlić powiększoną wersję tej ryciny.

figure-results-17
Rysunek 17Diagram oddziaływań z symulacji dynamiki molekularnej wolnego białka Etravirine oraz enumerowanego Etravirine. (A) Schemat oddziaływań dynamiki molekularnej wolnego białka. (B) Diagram oddziaływań dynamiki molekularnej Etravirine. (C) Schemat oddziaływań dynamiki molekularnej dla wyliczonej etrawiryny. Kliknij tutaj, aby wyświetlić powiększoną wersję tej ryciny.

figure-results-18
Rycina 18: Histogram kontaktów interakcyjnych między etrawiryną a białkiem. (A) Oś czasu kontaktów białko-ligand dla etrawiryny. (B) Oś czasu interakcji białko-ligand w funkcji czasu, obejmująca wiązania wodorowe, kontakty hydrofobowe, jonowe oraz mostki wodne. (C) Schemat przedstawiający szczegółowe interakcje między atomami ligandu a resztami białka. Wyświetlono jedynie interakcje występujące przez więcej niż 30% czasu symulacji (0,00 do 200 ns). Aby zobaczyć powiększoną wersję tej ryciny, kliknij tutaj.

figure-results-19
Rysunek 19: Histogram kontaktów interakcyjnych pomiędzy ponumerowaną etrawiryną a białkiem. (A) Oś czasu kontaktów białko-ligand dla ponumerowanej etrawiryny. (B) Oś czasu interakcji białko-ligand w czasie, obejmująca wiązania wodorowe, kontakty hydrofobowe, jonowe oraz mostki wodne. (C) Schemat przedstawiający szczegółowe interakcje pomiędzy atomami ligandu a resztami białka. Wyświetlono jedynie interakcje występujące przez ponad 30% czasu symulacji (0,00 do 200 ns). Kliknij tutaj, aby zobaczyć powiększoną wersję tego rysunku.

KlasaBiałko doceloweZakres związkówCałkowita liczba wybranych związków
1NL4-3 wild-type HIV-1(8a1-8e5 – EC50 (nM)a)25
2IIIB WT HIV-113a1-13d6 - EC50 (nM)a23
3RES056 szczep oporny na NNRTI13a1-13d6 - EC50 (nM)a23
4ROD szczep HIV-213a1-13d6 - EC50 (nM)a23
Suma związków wybranych do modelowania QSAR94

Tabela 1: Podsumowanie kryteriów klasyfikacji związków dla modelowania QSAR. Zbiór danych obejmujący 94 związki dihydrofuro[3,4-d] pirymidyny został zsyntetyzowany przez Kanga i współpracowników50 w celu ukierunkowania na różne szczepy HIV, mianowicie: HIV-1 typu dzikiego NL4-3, HIV-1 WT IIIB, szczep ROD HIV-2 oraz szczep RES056 oporny na NNRTI. Pochodne pirymidyny te pogrupowano w cztery klasy w zależności od białka docelowego w ramach przygotowań do przeprowadzenia szkolenia QSAR. Tabela podsumowuje sposób pogrupowania tych cząsteczek w cztery klasy na podstawie ich eksperymentalnych wartości EC50, które reprezentują potęgę działania leku, zdefiniowaną jako stężenie, przy którym lek wykazuje 50% swojego maksymalnego efektu. Skróty: NNRTI = nienukleozydowe inhibitory odwrotnej transkryptazy; QSAR = ilościowa zależność struktura-aktywność.

Nazwa białkaLigandWynik dokowania
1HQUEfavirenz-10.432
Etrawiryna-9.647
HBY_561-9.242
Dorawiryna-9.04
Newirapina-8.825
Rilpiwiryna-7.722
Delawirdyna-6.519

Tabela 2: Wyniki dokowania sześciu zoptymalizowanych NNRTIs oraz białka HIV-1.Wyniki dokowania dotyczą sześciu badanych NNRTIs. Bardziej ujemny wynik wskazuje na lepszą skuteczność wiązania liganda z białkiem. Efavirenz i Etravirine wykazały najkorzystniejsze wyniki, wynoszące odpowiednio -10,432 eV i -9,647 eV w stosunku do wyniku dokowania współkrystalizowanego liganda HBY561 (-9,242). Potencjalnymi ligandami były te, których wynik dokowania był bardziej ujemny niż -9,242 eV.

Klasa lekuDopasowanie aktywności przy 85 procentachKolumna1Kolumna2Kolumna3Kolumna4Kolumna5
WynikSDR2RMSEQ2Q2 MW (Hipoteza zerowa)
10.82230.32680.8150.24790.81850.1462
20.56710.29960.50670.19960.22640.4736
30.81720.36920.81140.15360.9065-1.1532
40.66730.3670.64360.17090.88520.1445
*SD – odchylenie standardowe,
 R2 – korelacja w zbiorze treningowym między rzeczywistymi a przewidywanymi wartościami aktywności, 
Q2 – korelacja między rzeczywistą a przewidywaną aktywnością w zbiorze testowym.
RMSE – pierwiastek błędu średniokwadratowego

Tabela 3: Parametry statystyczne modelu 2D-QSAR. Tabela przedstawia odchylenie standardowe, korelację w zbiorze treningowym między rzeczywistymi a przewidywanymi wartościami aktywności (R2) oraz wysoką wartość korelacji między rzeczywistą a przewidywaną aktywnością w zbiorze testowym dla każdej klasy (Q2). Wyższa wartość (R2) reprezentuje daną klasę.

LIGANDHOMOLUMOEprzerwa
Efavirenz-0.2242-0.069430.15477
Etrawiryna-0.21408-0.081410.13267
Newirapina-0.20602-0.077590.12843
HBY_561-0.19687-0.07480.12207
Delawirdyna-0.19651-0.079580.11693
Dorawiryna-0.22413-0.109160.11497
Rilpiwiryna-0.21343-0.112160.10127

Tabela 4: Przerwa energetyczna HOMO-LUMO dla zoptymalizowanych NNRTIs.Tabela przedstawia wyniki przerw energetycznych HOMO-LUMO uzyskane po optymalizacji sześciu NNRTIs.

LigandZoptymalizowane ligandyWymienione ligandy
Dorawiryna7.2297.374
Rilpiwiryna7.3027.279
Etrawiryna7.2297.374
Efavirenz7.2297.323
Delawirdyna7.3027.302
Newirapina7.2296.988
HBY_5617.2297.323

Tabela 5: Przewidywane wyniki aktywności zoptymalizowanych NNRTIs w porównaniu do odpowiadających im wyliczonych NNRTIs. Tabela porównuje przewidywane wyniki aktywności pomiędzy ligandami zoptymalizowanymi a wyliczonymi. Im wyższe wyniki aktywności, tym lepszy związek jako potencjalny kandydat na lek.

LIGANDPrzerwa energetyczna po optymalizacjiPrzerwa energetyczna po enumeracjiRóżnica luk między związkami zoptymalizowanymi a wyliczonymi
Dorawiryna0.1150.1260.011
Rilpiwiryna0.1010.1270.030
Etrawiryna0.1330.1260.010
Efavirenz0.1550.1260.030
Delawirdyna0.1170.1260.010
Newirapina0.1280.1260.002
HBY_5610.1220.1260.004

Tabela 6: Porównanie przewidywanych przerw energetycznych zenumerowanych NNRTIs i pierwotnej przerwy energetycznej zoptymalizowanych NNRTIs. Tabela przedstawia porównanie przerwy energetycznej HOMO-LUMO między zoptymalizowanymi a zenumerowanymi ligandami oraz różnice w przerwach energetycznych między nimi.

Wykatalogowany ligandWłaściwości reguły 5Kolumna 1Kolumna 2Kolumna 3Kolumna 4Kolumna 5Dodana grupa bocznaWymieniona ocena dokowaniaoryginalna wartość dokowaniaPonownie dokowane ligandy enumerowane
AlogPPSAHBDHBAmasa cząsteczkowaMPO
Dorawiryna2.4125.928441.80.49Hydroksylowy-8.894-9.04-9.739
Etrawiryna4.4140.938451.30.37Hydroksylowy-10.258-9.647-10.517
Efavirenz3.764.323330.70.75Amina-10.284-10.432-11.025
Newirapina2.658.114284.30.73Fluorek-9.112-8.825-9.445
HBY_5612.493.926358.50.66Amid-9.596-9.242-10.1
AlogP – (obliczony logarytm współczynnika podziału oktanol-woda).
PSA – (polarna powierzchnia cząsteczki)
HBD – (donory wiązań wodorowych)
HBA – (akceptory wiązań wodorowych)
MW – (masa cząsteczkowa)
MPO – (optymalizacja wieloparametrowa)

Tabela 7: Porównanie wyników dokowania pomiędzy oryginalnymi a wyliczonymi NNRTIs. Tabela przedstawia porównanie wyliczonych wyników dokowania, które są wynikami uzyskanymi po enumeracji. Oryginalne wyniki dokowania to wartości dla zoptymalizowanych NNRTIs. Wyniki redokowania wyliczonych ligandów to wartości dla ligandów, które zostały poddane enumeracji. Ponieważ proces enumeracji generuje predykcyjny wynik dokowania, ligandy musiały zostać poddane ponownemu dokowaniu tą samą metodą, co ligandy zoptymalizowane. Grupy dodane to te, które zostały przyłączone do ligandów podczas procesu enumeracji.

LigandΔGbindΔGCoulombΔGCovalentΔGHbondΔGLipoΔGpackΔGSolvΔGVdW
Etrawiryna-80.551-11.1962.491-1.509-27.447-4.82626.605-64.669
Etrawiryna enumerowana-89.684-17.9762.807-2.541-27.652-4.21126.034-66.146
HBY561-79.664-12.2610.994-0.521-26.052-1.27818.624-59.169
HBY561 enumerowany-82.719-13.4431.534-0.603-27.053-1.21020.600-62.544
Efavirenz-71.372-12.9841.151-0.834-25.353-1.37914.472-46.44
Efavirenz enumerowany-79.125-18.6021.753-2.159-25.409-1.19816.619-50.129

Tabela 8: Wyniki MMGBSA reXts wybranych ligandów dla 1HQU. Tabela przedstawia obliczenia mechaniki molekularnej z uogólnionym modelem Borna i polem powierzchni, które wykazują średnią wolną energię wiązania (Δ$\Delta G_{\text{bind}}$ kompleksu białko-ligand. Tabela przedstawia porównanie pierwotnie zoptymalizowanych ligandów, które wykazują dobre wyniki dokowania w stosunku do ligandów krystalicznych, oraz ich odpowiedników z enumeracji. Wybrany związek o wyższym wyniku dokowania i niskich fluktuacjach RMSD w symulacji dynamiki molekularnej podczas etapu ekwilibracji powinien również spełniać wymogi obliczeń MMGBSA, wykazując najbardziej ujemnąwolna energia wiązania (ΔGwiązać). W tym przypadku wyliczony etrawiryn spełnia oba parametry obliczeniowe.

Plik uzupełniający 1: Inne rXts uzyskane w niniejszym badaniu. Kliknij tutaj, aby pobrać ten plik.

Dyskusja

Metodę DFT26 wykorzystano do analizy właściwości elektronowych i stabilności różnych pochodnych pirymidyny, co nie ma miejsca w podobnych badaniach. Zastosowanie DFT pozwala na głębsze zrozumienie oddziaływań molekularnych na poziomie kwantowym, usprawniając proces projektowania NNRTI. W tym badaniu wybraliśmy sześć NNRTI o nieznanych wartościach aktywności. Pobraliśmy ich struktury molekularne z bazy danych PubChem i przeprowadziliśmy optymalizację geometryczną przy użyciu funkcjonalnego zestawu Minnesota 15 Local (MN15-L)32 i zestawu podstawowego 6-31++G (d,p)2 w pakiecie oprogramowania mechaniki kwantowej Gaussian 16 revision C0133 . Aby skonfigurować pliki wejściowe symulacji, użyliśmy GaussView 6.026. Struktury kwantowe zostały następnie zadokowane w miejscu aktywnym białka HIV-1.

W badaniu wykorzystano podejście QSAR51, które jest dostosowane do przewidywania aktywności biologicznej NNRTI na podstawie ich struktury chemicznej. Model ten został wygenerowany przy użyciu zestawu danych 94 pochodnych dihydrofuro[3,4-d] pirymidyny, co umożliwia identyfikację kluczowych grup funkcyjnych, które wpływają na aktywność przeciwwirusową. Włączenie solidnego modelu QSAR jest kluczowym postępem, ponieważ pomaga filtrować związki na wczesnym etapie procesu projektowania, zwiększając skuteczność odkrywania leków.

Wzięto pod uwagę przewidywaną aktywność nowych związków, a wyniki dokowania zostały zweryfikowane. Zaawansowana technika dokowania molekularnego pozwoliła na zbadanie interakcji wiążących między NNRTI a enzymem odwrotnej transkryptazy HIV-1. Uzupełnieniem były symulacje dynamiki molekularnej prowadzone przez dłuższy okres czasu (200 ns), które dostarczyły informacji na temat stabilności i zachowania kompleksów lek-białko w warunkach fizjologicznych. W ramach projektu Molecular Dynamics52 potwierdzono wyniki badań dokowania. Symulując zachowanie kompleksów lek-enzym w czasie, możemy ocenić dynamiczną stabilność i interakcje, które mogą nie zostać uchwycone w badaniach statycznego dokowania53. Podejście to zapewnia bardziej realistyczną ocenę tego, jak NNRTI mogą działać w systemach biologicznych. Symulacje wykazały, że etrawiryna powodowała wahania RMSD o około 4,5 A, podczas gdy wyliczana etrawiryna powodowała wahania RMSD o 3,5 A. Stwierdzono różne podobieństwa między etrawiryną a etryminowaną etrawiryną, przy czym wymieniony związek posiadał wzmocnione interakcje aminokwasowe w miejscu aktywnym białka. Niższy RMSD, ulepszone interakcje aminokwasowe i najwyższa wiązanie wolnej energii sugerują, że wymieniona etrawiryna może służyć jako realna alternatywa w leczeniu HIV/AIDS.

Na uwagę zasługuje wykorzystanie technik enumeracyjnych54 do generowania stereoizomerów i wykonywania skanowania izosterowego. Metoda ta pozwala naukowcom na zbadanie szerszej przestrzeni chemicznej dla potencjalnych NNRTI poprzez systematyczne modyfikowanie istniejących struktur, co może prowadzić do odkrycia związków o zwiększonej skuteczności i zmniejszonej oporności.

Badania łączą analizę HOMO-LUMO pochodzącą z obliczeń mechaniki kwantowej w celu oceny reaktywności i stabilności związków. Ten aspekt jest często pomijany w podobnychbadaniach55, a jednak dostarcza cennych informacji na temat właściwości elektronicznych, które mogą wpływać na aktywność biologiczną

Zastosowanie metody MMGBSA w tym badaniu do oszacowania wolnych energii wiązania wyliczonych NNRTI jako potencjalnych inhibitorów odwrotnej transkryptazy HIV-1 (RT) przedstawia nowy aspekt, który zwiększa znaczenie i potencjalny wpływ tej pracy na rozwój leczenia HIV.

Zastosowanie MMGBSA w tym badaniu zapewnia dokładniejsze oszacowanie wolnych energii wiązania dla wyliczonej etrawiryny w porównaniu z jej odpowiednikiem. Obliczając wartości wolnej energii wiązania, w badaniu ilościowo oceniano powinowactwo wiązania etrygowanej etrawiryny, które było znacznie wyższe (9,133 kcal/mol) niż w przypadku standardowej etrawiryny. Ten poziom szczegółowości wiążących obliczeń energetycznych pozwala na bardziej zniuansowane zrozumienie, w jaki sposób modyfikacje strukturalne mogą wpływać na skuteczność NNRTI, czego często brakuje w podobnych badaniach, które mogą opierać się wyłącznie na ocenach jakościowych.

Podana wartośćkowalencyjna ΔG wynosząca 1,477 kcal/mol dodatkowo wskazuje na siłę oddziaływania dla wymienionego związku. Podejście porównawcze jest innowacyjne, ponieważ nie tylko identyfikuje obiecującego kandydata, ale także sytuuje go w szerszym kontekście istniejących terapii, stanowiąc punkt odniesienia dla przyszłego rozwoju NNRTI.

Porównanie przedstawione w Tabeli 8 wskazuje również, że etrywata etrawiryna może być potencjalnym związkiem do zastosowania jako alternatywna ART, ponieważ jej swobodna energia wiązania (wiązanie ΔG) jest najwyższa w porównaniu z etrawiryną, HBY56, etrawiryną i ich wyliczonymi odpowiednikami.

Wyniki naszego badania mogą znacząco wpłynąć na opracowanie nowych środków terapeutycznych w leczeniu zakażenia wirusem HIV. Dodatkowo, takie podejście może pozwolić na zidentyfikowanie najbardziej obiecującego kandydata na anty-HIV. Odkrycia te mogą utorować drogę do racjonalnego projektowania nowych NNRTI HIV-1.

Włączenie MMGBSA w połączeniu z innymi metodami obliczeniowymi (takimi jak dokowanie molekularne, modelowanie QSAR i symulacje dynamiki molekularnej) zwiększa solidność wyników. To zintegrowane podejście pozwala na kompleksową ocenę związków, od optymalizacji strukturalnej po dynamiczne zachowanie w kontekście biologicznym. Taka metodologia jest stosunkowo nowa w badaniach NNRTI, gdzie badania często koncentrują się na wyizolowanych metodach bez holistycznego spojrzenia na proces projektowania leków. Badania wskazują na kilka obiecujących wyników i lepsze zrozumienie interakcji molekularnych między NNRTI a enzymem odwrotnej transkryptazy, co może stanowić podstawę przyszłych prac nad projektowaniem leków.

Podczas gdy inne badania w dziedzinie rozwoju NNRTI często koncentrują się na pojedynczych metodach obliczeniowych lub podstawowych analizach dokowania, badania te wyróżniają się połączeniem obliczeń chemii kwantowej z dynamiką molekularną i modelowaniem QSAR, zastosowaniem systematycznego podejścia do wyliczania, które rozszerza potencjalną przestrzeń chemiczną dla NNRTI, oraz wykorzystaniem kompleksowego in-silico ramy, które nie tylko przewidują powinowactwo wiązania, ale także oceniają stabilność i dynamikę interakcji lek-enzym.

Ogólnie rzecz biorąc, nowatorstwo tych badań polega na ich wieloaspektowym podejściu, wykorzystującym zaawansowane techniki obliczeniowe w celu sprostania wyzwaniom związanym z lekoopornością w leczeniu HIV, torując w ten sposób drogę do opracowania skuteczniejszych NNRTI. Odkrycia podkreślają potencjał opracowania skuteczniejszych NNRTI, które mogą przezwyciężyć ograniczenia obecnych terapii.

Niektóre przyszłe zastosowania podejścia opisanego w tym protokole obejmują następujące elementy.

Integracja z uczeniem maszynowym

Przyszłe badania mogą zintegrować MMGBSA z algorytmami uczenia maszynowego w celu poprawy mocy predykcyjnej wiązania szacunków energii swobodnej. Trenując modele na dużych zbiorach danych, naukowcy mogą zwiększyć dokładność i wiarygodność prognoz, umożliwiając lepszą identyfikację obiecujących kandydatów na NNRTI.

Zastosowanie w projektowaniu leków opartych na fragmentach

MMGBSA może być skutecznie wykorzystywana w projektowaniu leków opartych na fragmentach, gdzie małe fragmenty molekularne są optymalizowane w celu poprawy powinowactwa wiązania. Takie podejście może doprowadzić do identyfikacji nowych NNRTI o zwiększonej sile działania i selektywności wobec HIV-1.

Badanie mechanizmów oporności na leki

Technika ta może być zastosowana do badania interakcji wiązania NNRTI ze zmutowanymi szczepami HIV-1. Porównując wolne energie wiązania NNRTI zarówno z wariantami typu dzikiego, jak i opornymi, naukowcy mogą uzyskać wgląd w mechanizmy oporności na leki i pomóc w projektowaniu inhibitorów nowej generacji.

Ocena właściwości farmakokinetycznych

MMGBSA może być również stosowana do oceny właściwości farmakokinetycznych NNRTI, takich jak rozpuszczalność i przepuszczalność. Dzięki zrozumieniu, w jaki sposób te właściwości korelują z powinowactwem wiązania, naukowcy mogą zoptymalizować kandydatów na leki w celu uzyskania lepszych profili terapeutycznych.

Oświadczenia

Autorzy oświadczają, że nie są im znane żadne konkurencyjne interesy finansowe ani powiązania osobiste, które mogłyby mieć wpływ na pracę opisaną w tym artykule.

Podziękowania

Autorzy chcieliby podziękować Centrum Obliczeń o Wysokiej Wydajności (CHPC) za udostępnienie zasobów obliczeniowych oraz Wydziałowi Nauk Chemicznych Uniwersytetu w Johannesburgu.

Materiały

Lista materiałów użytych w tym artykule
NazwaFirmaNumer katalogowyKomentarze
GaussViewGaussViewV6.1.1
KNIME KNIMEV4.7.1
Schrödinger Maestro V13.6SCHRODINGER INC.Rok wydania licencji 2023-2

Bibliografia

  1. Fauci, A. S., Lane, H. C. Four decades of HIV/AIDS - much accomplished, much to do. N Engl J Med. 383 (1), 1-4 (2020).
  2. Kumar, V., Kishor, S., Ramaniah, L. M. Chemical reactivity analysis of deoxyribonucleosides and deoxyribonucleoside analogues (NRTIs): A first-principles density functional approach. J Mol Model. 18, 3969-3980 (2012).
  3. Shafer, R. W., Vuitton, D. A. Highly active antiretroviral therapy (haart) for the treatment of infection with human immunodeficiency virus type 1. Biomed Pharmacother. 53 (2), 73-86 (1999).
  4. Ginat, D. T., Schaefer, P. W. Highly active antiretroviral therapy (HAART). Neuroimaging Pharmacopoeia. Ginat, D. T., Small, D. T., Schaefer, P. W. , 229-238 (2022).
  5. Murray, C. J. L., et al. Global burden of bacterial antimicrobial resistance in 2019: A systematic analysis. Lancet. 399 (10325), 629-655 (2022).
  6. Razzaque, M. S. Commentary: Microbial resistance movements: An overview of global public health threats posed by antimicrobial resistance, and how best to counter. Front Public Health. 8, 629120-629120 (2021).
  7. Global HIV & AIDS statistics - Fact sheet. , UNAIDS. https://www.unaids.org/en/resources/fact-sheet (2023).
  8. Li, D., et al. HIV-1 pretreatment drug resistance and genetic transmission network in the southwest border region of China. BMC Infect Dis. 22 (1), 741(2022).
  9. Fact sheet: HIV drug resistance. , WHO. https://www.who.int/news-room/fact-sheets/detail/hiv-drug-resistance (2023).
  10. Hunt, G. M., et al. Prevalence of HIV-1 drug resistance amongst newly diagnosed HIV-infected infants age 4-8 weeks, enrolled in three nationally representative PMTCT effectiveness surveys, South Africa: 2010, 2011-12 and 2012-13. BMC Infect Dis. 19 (Suppl 1), 787-787 (2019).
  11. Koay, W. L. A., Kose-Otieno, J., Rakhmanina, N. HIV drug resistance in children and adolescents: Always a challenge. Curr Epidemiol Rep. 8 (3), 97-107 (2021).
  12. Wang, Z., et al. Contemporary medicinal chemistry strategies for the discovery and development of novel HIV-1 non-nucleoside reverse transcriptase inhibitors. J Med Chem. 65 (5), 3729-3757 (2022).
  13. Mbayo, V., Sookan, T. Correction to: Effects of a resistance training programme in people living with HIV in Zimbabwe. Sport Sci Health. 16 (4), 775(2020).
  14. Krawczyk, C. S., et al. Factors associated with delayed initiation of HIV medical care among infected persons attending a Southern HIV/AIDS clinic. South Med J. 99 (5), 472-481 (2006).
  15. Namasivayam, V., et al. The journey of HIV-1 non-nucleoside reverse transcriptase inhibitors (NNRTIs) from lab to clinic. J Med Chem. 62 (10), 4851-4883 (2019).
  16. Patel, P. H., Zulfiqar, H. Reverse transcriptase inhibitors. StatPearls. , StatPearls Publishing. Treasure Island (FL). (2023).
  17. Tan, J. J., et al. Therapeutic strategies underpinning the development of novel techniques for the treatment of HIV infection. Drug Discov Today. 15 (5-6), 186-197 (2010).
  18. Liu, N., et al. Novel HIV-1 non-nucleoside reverse transcriptase inhibitor agents: Optimization of diarylanilines with high potency against wild-type and rilpivirine-resistant e138k mutant virus. J Med Chem. 59 (8), 3689-3704 (2016).
  19. Anta, L., et al. Rilpivirine resistance mutations in HIV patients failing non-nucleoside reverse transcriptase inhibitor-based therapies. AIDS. 27 (1), 81-85 (2013).
  20. Sarafianos, S. G., et al. Structure and function of HIV-1 reverse transcriptase: Molecular mechanisms of polymerization and inhibition. J Mol Biol. 385 (3), 693-713 (2009).
  21. Vingerhoets, J., et al. Tmc125 displays a high genetic barrier to the development of resistance: Evidence from in vitro selection experiments. J Virol. 79 (20), 12773-12782 (2005).
  22. Ripamonti, D., Bombana, E., Rizzi, M. Rilpivirine: Drug profile of a second-generation non-nucleoside reverse transcriptase hiv-inhibitor. Expert Rev Anti-infect Ther. 12 (1), 13-29 (2014).
  23. Lambert-Niclot, S., et al. Prevalence of pre-existing resistance-associated mutations to rilpivirine, emtricitabine and tenofovir in antiretroviral-naive patients infected with B and non-B subtype HIV-1 viruses. J Antimicrob Chemother. 68 (6), 1237-1242 (2013).
  24. Wainberg, M. A. Combination therapies, effectiveness, and adherence in patients with HIV infection: Clinical utility of a single tablet of emtricitabine, rilpivirine, and tenofovir. HIV AIDS (Auckl). 5, 41-49 (2013).
  25. Namasivayam, V., et al. The journey of HIV-1 non-nucleoside reverse transcriptase inhibitors (NNRTIs) from lab to clinic. J Med Chem. 62 (10), 4851-4883 (2018).
  26. Jordaan, M. A., Ebenezer, O., Damoyi, N., Shapi, M. Virtual screening, molecular docking studies and DFT calculations of FDA approved compounds similar to the non-nucleoside reverse transcriptase inhibitor (NNRTI) efavirenz. Heliyon. 6 (8), e04642(2020).
  27. Soltani, A., et al. Application of molecular docking for the development of improved HIV-1 reverse transcriptase inhibitors. Curr Comput Aided Drug Des. 17 (4), 538-549 (2021).
  28. Vanangamudi, M., Palaniappan, S., Kathiravan, M. K., Namasivayam, V. Strategies in the design and development of non-nucleoside reverse transcriptase inhibitors (NNRTIs). J Viruses. 15 (10), 1992(2023).
  29. Mohapatra, R. K., et al. Computational investigations of three main drugs and their comparison with synthesized compounds as potent inhibitors of SARS-CoV-2 main protease (Mpro): DFT, QSAR, molecular docking, and in silico toxicity analysis. J King Saud Univ Sci. 33 (2), 101315-101315 (2021).
  30. Kim, S., et al. Pubchem substance and compound databases. Nucleic Acids Res. 44 (D1), D1202-D1213 (2016).
  31. Ofem, M. I., et al. Synthesis, spectral characterization, and theoretical investigation of the photovoltaic properties of (E)-6-(4-(dimethylamino) phenyl) diazenyl)-2-octyl-benzoisoquinoline-1, 3-dione. BMC Chem. 16 (1), 109(2022).
  32. Yu, H. S., He, X., Truhlar, D. G. MN15-L: A new local exchange-correlation functional for kohn-sham density functional theory with broad accuracy for atoms, molecules, and solids. J Chem Theory Comput. 12 (3), 1280-1293 (2016).
  33. Sasitha, T., John, W. J. Design, docking, and DFT investigations of 2,6-bis(3,4-dihydroxyphenyl)-3-phenethylpiperidin-4-one. Heliyon. 7 (2), e06127-e06127 (2021).
  34. Maurya, S. K., Maurya, A. K., Mishra, N., Siddique, H. R. Virtual screening, ADME/T, and binding free energy analysis of anti-viral, anti-protease, and anti-infectious compounds against nsp10/nsp16 methyltransferase and main protease of SARS CoV-2. J Recept Signal Transduct Res. 40 (6), 605-612 (2020).
  35. Singh, A. K., et al. Current insights and molecular docking studies of HIV-1 reverse transcriptase inhibitors. Chem Biol Drug. 103 (1), e14372(2024).
  36. Kralj, S., Jukič, M., Bren, U. Comparative analyses of medicinal chemistry and cheminformatics filters with accessible implementation in konstanz information miner (KNIME). Int J Mol Sci. 23 (10), 5727(2022).
  37. Bastikar, V., Bastikar, A., Gupta, P. P. Quantitative structure-activity relationship-based computational approaches. Computational Approaches for Novel Therapeutic and Diagnostic Designing to Mitigate SARS-CoV-2 Infection. , 191-205 (2022).
  38. Mazouin, B., Schöpfer, A. A., Von Lilienfeld, O. A. Selected machine learning of HOMO-LUMO gaps with improved data-efficiency. Mater Adv. 3 (22), 8306-8316 (2022).
  39. Singh, V. K., et al. In silico design, synthesis and anti-HIV activity of quinoline derivatives as non-nucleoside reverse transcriptase inhibitors (NNRTIs). Comput Biol Chem. 98, 107675(2022).
  40. Lanka, G., et al. Pharmacophore-based virtual screening, 3D QSAR, docking, ADMET, and MD simulation studies: An in silico perspective for the identification of new potential HDAC3 inhibitors. Comput Biol Med. 166, 107481(2023).
  41. Singh, K. D., Muthusamy, K. Molecular modeling, quantum polarized ligand docking and structure-based 3D-QSAR analysis of the imidazole series as dual AT1 and ETA receptor antagonists. Acta Pharmacol Sin. 34 (12), 1592-1606 (2013).
  42. Salo-Ahen, O. M., et al. Molecular dynamics simulations in drug discovery and pharmaceutical development. Processes. 9 (1), 71(2020).
  43. Raniolo, S., Limongelli, V. Improving small-molecule force field parameters in ligand binding studies. Front Mol Biosci. 8, 760283(2021).
  44. Lu, C., et al. OPLS4: Improving force field accuracy on challenging regimes of chemical space. J Chem Theory Comput. 17 (7), 4291-4300 (2021).
  45. Goswami, N., Singh, A., Bharadwaj, S., Sahoo, A. K., Singh, I. K. Targeting neuroblastoma by small-molecule inhibitors of human ALYREF protein: Mechanistic insights using molecular dynamics simulations. J Biomol Struct Dyn. 42 (3), 1352-1367 (2023).
  46. Chirico, N., Gramatica, P. Real external predictivity of QSAR models. Part 2. New intercomparable thresholds for different validation criteria and the need for scatter plot inspection. J Chem Inf Model. 52 (8), 2044-2058 (2012).
  47. Tabti, K., Sbai, A., Maghat, H., Lakhlifi, T., Bouachrine, M. Computational exploration of the structural requirements of triazole derivatives as colchicine binding site inhibitors. ChemistrySelect. 8 (26), e202301707(2023).
  48. Miar, M., Shiroudi, A., Pourshamsian, K., Oliaey, A. R., Hatamjafari, F. Theoretical investigations on the HOMO-LUMO gap and global reactivity descriptor studies, natural bond orbital, and nucleus-independent chemical shifts analyses of 3-phenylbenzo[D]thiazole-2(3H)-imine and its para-substituted derivatives: Solvent and substituent effects. J Chem Res. 45 (1-2), 147-158 (2020).
  49. Saldivar-Gonzalez, F., Huerta-García, C., Medina-Franco, J. Chemoinformatics-based enumeration of chemical libraries: A tutorial. J Cheminf. 12, 64-64 (2020).
  50. Kang, D., et al. Identification of dihydrofuro[3,4-d]pyrimidine derivatives as novel HIV-1 non-nucleoside reverse transcriptase inhibitors with promising antiviral activities and desirable physicochemical properties. J Med Chem. 62 (3), 1484-1501 (2019).
  51. Viira, B., Garcia-Sosa, A. T., Maran, U. Chemical structure and correlation analysis of HIV-1 NNRT and NRT inhibitors and database-curated, published inhibition constants with chemical structure in diverse datasets. J Mol Graph. 76, 205-223 (2017).
  52. Javed, M. R. CADD and molecular dynamic simulations: Potential impacts to conventional medicines. Comb Chem High Throughput Screen. 25 (4), 658-659 (2022).
  53. Jiang, X., et al. Exploiting the tolerant region i of the non-nucleoside reverse transcriptase inhibitor (NNRTIi) binding pocket. Part 2: Discovery of diarylpyrimidine derivatives as potent HIV-1 NNRTIs with high FSP3 values and favorable drug-like properties. Eur J Med Chem. 213, 113051(2021).
  54. Zhang, T., Jiang, S., Li, T., Liu, Y., Zhang, Y. Identified isosteric replacements of ligands' glycosyl domain by data mining. ACS Omega. 8 (28), 25165-25184 (2023).
  55. Sule, L., Gupta, S., Jain, N., Sapre, N. S. In silico induction of missense mutation in NNRTI protein: Computational modelling and stability study of modelled proteins. J Math Chem. 62, 2776-2797 (2024).

Przedruki i uprawnienia

Tagi

Dokowanie molekularnesymulacje dynamiki molekularnejteoria funkcjonału gęstościkompleks białko-ligandwolna energia wiązaniaoporność HIV-1 na lekipochodne pirymidynyprzygotowanie ligandu
Ilościowa zależność struktura-aktywność, przewidywanie aktywności i dynamika molekularna nienukleotydowych inhibitorów odwrotnej transkryptazy
Skontaktuj się z nami
Poleć bibliotece

Prywatność

Warunki korzystania z usługi

Zasady

Prawa autorskie © 2026 MyJoVE Corporation. Wszelkie prawa zastrzeżone.