Abstrak
Resistensi adaptif nongenetik terhadap kemoterapi, yang didorong oleh penataan ulang transkripsi, muncul sebagai mekanisme penting dalam kelangsungan hidup tumor. Dalam studi ini, kami menggabungkan transkriptomik longitudinal dengan analisis pola temporal untuk menyelidiki mekanisme khusus pasien yang mendasari resistensi yang didapat pada kanker payudara. Biopsi tumor yang cocok (praperawatan, pascaperawatan, dan normal yang berdekatan) dikumpulkan dari pasien kanker payudara yang menerima kemoterapi neoadjuvan. Transkriptom dianalisis dengan klasifikasi pola gen longitudinal untuk melacak perubahan ekspresi gen khusus pasien yang terjadi selama perawatan. Temuan kami mengungkapkan bahwa gen yang terkait dengan resistensi sudah mengalami disregulasi pada tumor primer, yang menunjukkan adanya keadaan toleran obat yang sudah ada sebelumnya. Sementara setiap pasien menunjukkan penataan ulang gen terkait resistensi yang unik, perubahan ini menyatu menjadi sejumlah modul fungsional yang tidak diatur. Khususnya, pasien yang menerima perawatan yang sama menunjukkan penataan ulang gen dan jalur yang berbeda, yang mengungkapkan rute paralel dan individual menuju resistensi. Sebagai kesimpulan, kami mengusulkan agar sel tumor bertahan hidup dari kemoterapi dengan mempertahankan atau memperkuat keadaan toleran obat yang sudah ada sebelumnya yang menghindari aksi obat. Kami menyarankan bahwa “peta kemoresistom” yang diindividualisasikan dapat mengidentifikasi kerentanan kanker dan menginformasikan strategi terapi yang dipersonalisasi untuk mengatasi atau mencegah resistensi.

Singkatan
ADC
konjugat obat antibodi
ECM
matriks ekstraseluler
FFPE
formalin-fiksasi parafin-tertanam
SDM+
reseptor hormon positif
IRB
dewan peninjau kelembagaan
Anggota Parlemen
Miller–Payne
Bahasa Indonesia: NAT
pengobatan neoadjuvant
Bahasa Inggris
jaringan payudara normal
RCB
beban sisa kanker
TNBBS
kanker payudara triple negatif
VST
transformasi stabilisasi varians
1 Pendahuluan
Meskipun resistensi obat tumor sangat penting terhadap morbiditas dan mortalitas kanker, pemahaman kita tentang berbagai mekanisme molekuler yang terlibat dalam resistensi terbatas [ [ 1 ] ]. Respons aktual dari pasien individu tetap menjadi ‘kotak hitam.’ Ahli onkologi tidak dapat secara akurat memprediksi tumor mana, yang memenuhi syarat untuk kemoterapi, yang akan dihilangkan oleh obat-obatan, dan resistensi intrinsik terhadap kemoterapi tetap menjadi teka-teki yang substansial. Baru-baru ini, seleksi klonal genetik yang berlaku, sebagai mekanisme utama resistensi, ditantang oleh indikasi yang muncul dari mekanisme adaptif nongenetik dari resistensi obat [ [ 2 , 3 ] ]. Bukti kumulatif menunjukkan bahwa sebagai respons terhadap pengobatan, tumor mengadopsi keadaan toleran obat melalui pemrograman ulang transkripsi [ [ 4 – 6 ] ]. Untuk bertahan hidup dari kemoterapi, sel kanker menyusun ulang jalur molekuler, sehingga terhindar dari obat antiproliferatif. Namun, karena kanker adalah proses yang berkelanjutan, rute individu menuju resistensi pada setiap jenis kanker dan setiap pasien sulit dipahami.

Kanker payudara saat ini merupakan jenis kanker wanita yang paling umum dan penyebab kedua kematian akibat kanker pada wanita [ [ 7 ] ]. Upaya internasional utama telah membuat profil tumor primer dan menentukan arsitektur molekuler kanker payudara [ [ 8 – 12 ] ]. Sekarang jelas bahwa kanker payudara adalah penyakit heterogen, diklasifikasikan menjadi beberapa subtipe: reseptor hormon positif (HR+), HER2-positif dan kanker payudara triple-negatif (TNBC). Sementara terapi target dan imunoterapi baru sedang dikembangkan, kemoterapi tetap menjadi andalan pengobatan untuk pasien kanker payudara berisiko tinggi awal, dan tulang punggung konjugat antibodi-obat (ADC) yang baru dikembangkan [ [ 13 ] ]. Dengan demikian, mengungkap mekanisme resistensi kemoterapi merupakan kebutuhan krusial.

Untuk menjelaskan mekanisme resistensi, cara yang ideal adalah melakukan studi longitudinal, dengan melihat tumor pra- dan pasca-perawatan yang cocok. Pengaturan pra operasi neoadjuvant (NAT) menyediakan kesempatan unik untuk mempelajari efek kemoterapi pada tumor di kehidupan nyata melalui keberadaan inheren biopsi diagnostik pra-perawatan dan spesimen operasi pasca-perawatan [ [ 14 ] ]. NAT digunakan secara luas untuk perawatan pasien kanker payudara berisiko tinggi dan memberikan informasi prediktif dan prognostik jangka panjang dengan mengevaluasi kemosensitivitas tumor [ [ 15 ] ]. Pencarian biomarker yang memprediksi respons atau target molekuler yang menghindari resistensi memicu banyak penelitian yang memeriksa tumor residual pasca-NAT. Sementara penelitian awal membuat profil sampel tumor primer satu titik waktu [ [ 16 ] ] atau tumor pasca-NAT residual [ [ 17 ] ], penelitian selanjutnya membandingkan tumor pra- dan pasca-perawatan [ [ 18 – 22 ] ], yang menghasilkan penemuan gen dan jalur yang diubah oleh pengobatan. Namun, sebagian besar penelitian didasarkan pada ekspresi diferensial statistik dan tidak dapat memberikan pandangan resistensi yang dipersonalisasi. Penelitian yang lebih baru memeriksa perubahan perjalanan waktu setelah berbagai modalitas pengobatan neoadjuvant, menggunakan pengambilan sampel serial yang cocok dari tiga titik waktu: pra-NAT, jangka pendek dan jangka panjang pasca pengobatan [ [ 6 , 14 , 23 – 26 ] ]; meskipun demikian, mekanisme yang mendasari resistensi nongenetik dan tingkat keragaman fenotipik antara pasien setelah pengobatan masih belum jelas. Pemrograman ulang transkripsi sebagai respons terhadap kemoterapi dapat dimediasi oleh mekanisme stokastik yang mendorong heterogenitas nongenetik [ [ 27 ] ]. Masih belum jelas, apakah resistensi nongenetik melibatkan pemrograman ulang aktif gen-gen spesifik selama pengobatan atau, sebagai alternatif, mengintensifkan gen yang disiapkan secara transkripsi. Untuk menyelidiki hal ini dalam sampel yang relevan secara klinis, kami mengumpulkan sekelompok pasien kanker payudara neoadjuvant dan membuat profil transkriptom mereka sebelum dan sesudah pengobatan, dibandingkan dengan epitel normal yang berdekatan dari pasien yang sama. Dengan demikian, kami memeriksa modulasi yang diinduksi obat relatif terhadap keadaan deregulasi tumor primer vs payudara normal yang cocok. Untuk memahami jalur pribadi menuju resistensi, kami menerapkan algoritma analisis pola yang berorientasi pada pasien untuk menangani kumpulan data longitudinal pada tingkat gen. Kami sebelumnya mengembangkan pendekatan analisis pola ini untuk menyelidiki modulasi ekspresi miRNA setelah kekambuhan [ [ 28 ]], dan selanjutnya menerapkannya pada dataset NAT proteomik longitudinal pasien kanker payudara, mengungkapkan keterlibatan biosintesis prolin dalam resistensi [ [ 29 ] ].

Analisis kumpulan data longitudinal kami menggambarkan prinsip-prinsip penting dalam mekanisme resistensi nongenetik adaptif. Dengan mendekonstruksi modulasi pola dinamis untuk setiap gen per pasien, kami mendeteksi pengkabelan ulang beberapa jalur secara bersamaan, meskipun dengan cara yang unik untuk setiap pasien. Pemetaan kemoresistom pribadi dapat memberikan pemahaman mekanistik yang mendalam tentang resistensi terapi dengan implikasi klinis di masa mendatang untuk desain rencana perawatan yang rasional.

2 Metode
2.1 Perakitan kelompok
Kohort terdiri dari 29 pasien kanker payudara yang menjalani terapi neoadjuvant, dan tidak mencapai respons lengkap, sehingga menghasilkan respons parsial dan tumor residual yang tersedia. Hanya pasien dengan respons parsial yang disertakan, karena kami ingin fokus pada perubahan dinamis antara tumor primer dan tumor residual. Informasi klinikopatologis untuk kohort dirangkum dalam Tabel 1 dan dirinci dalam Tabel S1 . Untuk setiap pasien, kami mengumpulkan triplet tumor praperawatan (biopsi diagnostik), tumor residual pascaperawatan (diambil sampelnya saat operasi, ~6–8 bulan setelah biopsi praperawatan), dan epitel normal yang berdekatan (diambil sampelnya dari spesimen operasi pascaperawatan). Sampel dikumpulkan sebagai jaringan formalin-fixed paraffin-embedded (FFPE) dari arsip patologi Sheba Medical Center. Yang penting, pengumpulan jaringan dan waktu fiksasi distandarisasi di departemen patologi. Selain itu, kami sebelumnya menunjukkan bahwa profil transkriptom dari jaringan FFPE menggunakan protokol kami yang dioptimalkan mirip dengan transkriptom dari jaringan snap-frozen yang sama [ [ 30 ] ]. Sampel epitel payudara normal tambahan dikumpulkan dari spesimen reduksi payudara dari enam individu sehat. Secara keseluruhan, total 97 sampel disertakan dalam penelitian ini. Dua ahli patologi independen memeriksa slide untuk menentukan selularitas tumor, menandai daerah untuk makrodiseksi untuk memperkaya selularitas tumor dalam ekstraksi RNA. Ahli patologi menilai respons tumor dengan skor Miller–Payne (MP) [ [ 31 ] ], yang menilai reduksi selularitas tumor dengan membandingkan sampel biopsi praperawatan dengan spesimen bedah pascaperawatan, dan dengan RCB (kelas Beban Kanker Residual) [ [ 32 ] ]. Kelangsungan hidup bebas kekambuhan menurut skor Miller–Payne dilakukan dengan analisis Kaplan–Meier, menghitung nilai P log-rank . Penelitian ini disetujui oleh Dewan Peninjauan Institusional (IRB) Sheba Medical Center (nomor persetujuan #8736-11), sesuai dengan prinsip-prinsip Deklarasi Helsinki, dan dengan Peraturan Kesehatan Masyarakat. Protokol penelitian mengizinkan pengecualian penuh untuk formulir persetujuan untuk sampel yang dianonimkan. Oleh karena itu, semua sampel dianonimkan sebagaimana didefinisikan dalam protokol penelitian.

Tabel 1. Parameter klinikopatologi untuk seluruh kelompok

LN, kelenjar getah bening; MP, Miller–Payne; RCB, beban kanker residual; RFS, kelangsungan hidup bebas kekambuhan. Kemoterapi: A – adriamycin (doxorubicin), C – cyclophosphamide, T – taxol (paclitaxel), H – herceptin (sebagian besar pasien HER2+ diobati sebelum Pertuzumab dimasukkan dalam pengobatan standar).
a Satu pasien memiliki tumor bilateral (HR+ diambil sampelnya dan HER2+ tidak diambil sampelnya), oleh karena itu menerima herceptin.

2.2 Persiapan pustaka mRNAseq
Sampel tumor FFPE dibedah secara makro untuk memperkaya selularitas tumor. Potongan 5-μm dideparafinisasi pada suhu 90 °C selama 5 menit. Total RNA serta DNA diekstraksi menggunakan kit isolasi asam nukleat (AllPrep, Qiagen, Hilden, Jerman) sesuai dengan instruksi protokol. Konsentrasi RNA ditentukan dengan fluorometer Qubit (Thermofisher Scientific, Waltham, MA, AS). Profil transkriptom keseluruhan dari sampel yang diarsipkan dicapai dengan prosedur mRNAseq yang andal dan hemat biaya yang sebelumnya kami optimalkan untuk sampel FFPE [ [ 30 ] ]. Secara singkat, pustaka mRNA-seq dibangun dengan kit persiapan sampel RNA Truseq v2 (Illumina, San Diego, CA, AS). Ukuran pustaka dievaluasi oleh Tapestation (Agilent, Santa Clara, CA, AS). Dimer primer dihilangkan menggunakan 1× Agencourt RNAClean XP Beads. Pustaka sekuensing dibangun dengan TruSeq SBS Kit dan dimultipleks untuk menjalankan 6–10 sampel per jalur dengan panjang pembacaan 60 bp pada satu ujung dalam instrumen Illumina HiSeq V4. Kedalaman sekuensing rata-rata adalah ~27 juta pembacaan per sampel.

2.3 Penyelarasan dan pemetaan urutan
Adaptor Illumina dan peregangan PolyA/T dipangkas menggunakan cutadapt [ [ 33 ] ]. Bacaan yang lebih pendek dari 30 bp dibuang. Bacaan dipetakan ke genom manusia (GRCh38) menggunakan star v2.4.2a [ [ 34 ] ]. Secara keseluruhan, pemetaan bacaan ke genom tinggi, dengan rata-rata 87% bacaan dipetakan dan 61% bacaan yang dipetakan secara unik dihitung. Median bacaan yang dipetakan secara unik adalah 10,7 juta bacaan. Penghitungan dilanjutkan pada gen yang dianotasi dalam Ensembl rilis 88, menggunakan htseq-count [ [ 35 ] ]. Dari 97 sampel, empat sampel normal gagal dalam langkah persiapan pustaka karena kuantitas RNA yang rendah. Enam sampel tambahan gagal dalam kontrol kualitas untuk pemetaan dan tiga sampel dikeluarkan karena selularitas tumor yang rendah. Secara keseluruhan, 84 sampel dari 29 pasien disertakan dalam analisis umum. Akhirnya, 23 pasien memiliki sampel triplet lengkap yang lolos semua QC data dan digunakan untuk analisis pola (Gbr. S1 ). Sebuah alur kerja dibangun menggunakan snakemake [ [ 36 ] ]. Analisis data sebagian besar dilakukan dengan bahasa R [ [ 37 ] ] menggunakan paket yang ditentukan.

2.4 Evaluasi data mRNAseq
Rentang dinamis data mRNAseq dikorelasikan dengan nilai skor patologis pewarnaan imunohistokimia untuk: ER, PR, HER2, dan KI67 jika tersedia (Gbr. S2A,B ). Persentase sel positif Ki67 dinilai secara buta terhadap data ekspresi menggunakan algoritma analisis gambar Virtuoso otomatis (Ventana Medical Systems, Tucson, AZ, AS). Rentang dinamis data mRNAseq pada skor patologis tertinggi untuk reseptor jauh lebih tinggi daripada tingkat pewarnaan imunohistokimia (Gbr. S2C ). Untuk mengevaluasi kemurnian sampel tumor dari hitungan mRNAseq, kami menggunakan algoritma ESTIMATE [ [ 38 ] ] dan membandingkannya dengan kemurnian tumor, yang diperkirakan oleh dua ahli patologi independen setelah makrodiseksi. Ada kesesuaian yang baik antara evaluasi patologis dan kemurnian yang dihitung oleh algoritma ESTIMATE.

2.5 Analisis ekspresi diferensial
Untuk menganalisis ekspresi diferensial antara tipe sampel yang cocok (normal—praperawatan—pascaperawatan), kami menggunakan paket deseq 2 [ [ 39 ] ] dengan betaPrior, cooksCutoff, dan parameter Filtering independen yang ditetapkan ke False. Informasi pasien ditambahkan ke model deseq 2 untuk mengendalikan efeknya. Nilai P mentah disesuaikan untuk pengujian berganda menggunakan prosedur Benjamini dan Hochberg. Untuk memperkirakan jumlah gen yang diekspresikan secara diferensial antara setiap pasangan tipe sampel, kontras berikut ditetapkan: tumor pascaperawatan vs normal, tumor praperawatan vs normal, dan pascaperawatan vs praperawatan. Gen-gen tersebut difilter untuk mempertahankan gen dengan perubahan lipatan absolut di atas 2, nilai P yang disesuaikan di bawah atau sama dengan 0,05, dan jumlah sedikitnya 100 dalam sedikitnya satu sampel. Sekitar 6900 gen melewati ambang batas ini. Subtipe molekuler dihitung dari data mRNAseq oleh paket genefu [ [ 40 ] ] menggunakan algoritma subtipe SCMOD2 [ [ 41 ] ] (dirinci dalam Tabel S1 ).

2.6 Klasifikasi pola
Untuk menganalisis modifikasi gen di tiga titik waktu (normal, praperawatan, dan pascaperawatan), kami mendefinisikan semua delapan pola teoritis yang mungkin (Gbr. 2C ) seperti yang diterbitkan sebelumnya [ [ 28 ] ]. Untuk menetapkan pola untuk setiap gen, kami mengorelasikan perubahan nilai ekspresi, dihitung dengan deseq 2 Variance Stabilizing Transformation (VST), ke pola teoritis. VST menghitung rata-rata geometrik sampel, memperkirakan dispersi untuk setiap gen, menyesuaikan tren dispersi rata-rata, dan kemudian mengubah data untuk menstabilkan varians di seluruh rata-rata. Setiap gen ditugaskan ke pola teoritis (menurut korelasi maksimal). Gen dianggap terklasifikasi menjadi satu pola utama jika diklasifikasikan ke dalam satu pola tunggal pada setidaknya 50% pasien. Gen dianggap terklasifikasi ke dalam pola bimodal jika gen diklasifikasikan ke masing-masing dari dua pola pada setidaknya 30% pasien dan masing-masing dari dua pola diamati setidaknya pada lima pasien, berdasarkan distribusi frekuensi pola.

2.7 Mendefinisikan gen yang terkait dengan respons
Untuk gen dengan pola utama yang bertahan: klasifikasi P3 atau P3 + P7 atau P4 atau P4 + P8 pada lebih dari 50% pasien. Untuk gen dengan pola bimodal, kami membandingkan skor Miller–Payne dari pasien dengan pola “respons baik” (P1 atau P2) versus skor pasien dengan pola “respons buruk” (P3, P7 atau P4, P8). Uji Wilcoxon ranking sum digunakan untuk menghitung signifikansi hubungan pola dengan respons. Uji ini dihitung untuk setiap klaster gen secara terpisah: P1 vs P3/P7 dan P2 vs P4/P8 (Gbr. 3B ). Uji ini diulang hanya untuk pasien HR+.

2.8 Algoritma Pathifier
Kami menggunakan algoritma Pathifier [ [ 42 , 43 ] ] untuk menghitung skor deregulasi jalur (PDS) untuk setiap sampel pasien untuk daftar gabungan basis data jalur: proses GO_biological (2018), KEGG (2019), dan Reactome (2016), semuanya diunduh dari situs web EnricheR [ [ 44 ] ]. Secara singkat, Pathifier menganalisis setiap jalur, dan menetapkan untuk setiap sampel i dan jalur P skor yang memperkirakan sejauh mana perilaku jalur P menyimpang dari sampel normal. Untuk menentukan PDS, kami menggunakan hitungan normal dari semua 97 sampel. Untuk setiap jalur, PDS didefinisikan sebagai jarak setiap sampel dari proyeksinya pada kurva utama yang dihitung dari semua sampel ke proyeksi sampel normal (Gbr. 1C ). Sampel normal didefinisikan sebagai sampel normal yang sehat dan sampel normal yang berdekatan.

Gbr. 1
Buka di penampil gambar
Kekuatan Gambar
Transkriptomik longitudinal dari sampel yang cocok dari pasien kanker payudara yang menerima kemoterapi neoadjuvan. (A) Sampel yang diarsipkan dari 29 pasien kanker payudara yang menjalani terapi neoadjuvan dikumpulkan. Dari setiap pasien kami menganalisis praperawatan, pascaperawatan, dan epitel payudara normal yang berdekatan yang cocok. Selain itu, sampel payudara normal dari tiga wanita sehat disertakan. (B) Analisis komponen utama dari 1000 gen paling bervariasi dari seluruh data mRNAseq transkriptom. Ekspresi gen dikelompokkan menurut jenis sampel. Sampel payudara normal yang berdekatan mirip dengan mamae normal dari wanita sehat. (C) Contoh penghitungan skor deregulasi jalur (PDS) menggunakan algoritma Pathifier untuk jalur yang dipilih dari kumpulan data KEGG. Kurva utama (hitam) melewati awan titik yang mewakili berbagai sampel. Sampel (diwarnai menurut jenisnya) diproyeksikan ke kurva. (D) Skor Deregulasi Jalur (PDS) dihitung untuk setiap sampel per jalur menggunakan algoritma Pathifier. Setiap baris mewakili jalur yang berbeda (KEGG, Reactome, dan GO), yang mana PDS relatif terhadap semua sampel normal (sehat + berdekatan) dihitung. Setiap kolom mewakili sampel yang berbeda pada masing-masing pasien. Khususnya, jaringan payudara normal yang sehat (NBR) sangat mirip dengan jaringan normal yang berdekatan (N), seperti yang ditunjukkan oleh PDS rendah secara keseluruhan (biru). Sampel tumor sebelum perawatan menunjukkan skor yang paling tidak teratur, seperti yang ditunjukkan oleh tingkat merah secara keseluruhan. Sampel tumor setelah perawatan bervariasi. (E) Modulasi temporal dalam skor deregulasi per pasien untuk dua jalur representatif. RCB, beban kanker residual.

2.9 Analisis pengayaan jalur dan perhitungan skor resistensi
Analisis pengayaan jalur untuk gen yang secara signifikan terkait dengan respons (uji Wilcoxon rank-sum, nilai- P < 0,05) dihitung menggunakan paket clusterprofiler [ [ 45 ] ], menggunakan jalur dari basis data berikut: proses GO_biological (2018), KEGG (2019), dan Reactome (2016) (diunduh dari situs web basis data EnrichR [ [ 44 ] ]). Untuk setiap jalur yang diperkaya, kami menghitung skor resistensi per pasien dan nilai- P . Skor resistensi didefinisikan sebagai fraksi “gen resistensi” (gen dengan pola P3, P4, P7, atau P8) dalam jalur tersebut. Hanya gen yang diekspresikan secara berbeda secara signifikan yang dipertimbangkan. Nilai- P untuk skor ini dihitung dengan mem-bootstrapping set acak gen dengan ukuran jalur yang sama (diambil sampelnya 10.000 kali). Untuk setiap pasien di setiap jalur, nilai- P dihitung dengan membagi jumlah skor acak, yang lebih besar atau sama dengan skor yang diamati, dengan jumlah permutasi. Untuk mengurangi bias subtipe, pengayaan jalur dilakukan dua kali: untuk seluruh kelompok pasien (23 pasien) dan untuk pasien HR+ (20 pasien), tidak termasuk tiga TNBC. Kami kemudian memotong jalur yang diperkaya antara kedua kelompok untuk mencapai daftar jalur yang kuat (33 jalur). Untuk mengimbangi keragaman pasien, kami menambahkan ke daftar bersama lima jalur paling signifikan yang unik untuk setiap kelompok, menghasilkan 43 jalur. Jalur yang diperkaya kemudian dikurasi ke daftar akhir 24 jalur, memperhitungkan signifikansi dan redundansi fungsional. Jalur ini dibagi menjadi modul fungsional. Untuk representasi sederhana modul resistensi pada setiap pasien, kami membuat peta kemoresistom pribadi. Barplot melingkar yang menggambarkan skor resistensi dan nilai- P di setiap modul diplot untuk setiap pasien.

2.10 Analisis STRING
Untuk mengidentifikasi jaringan fungsional gen resistensi (254 gen, nilai- P dengan uji peringkat Wilcoxon < 0,05), kami menjalankan analisis STRING (v. 11,5) [ [ 46 ] ]. Skor interaksi ditetapkan pada keyakinan tertinggi (0,9). Pengelompokan K -means (didefinisikan untuk 11 klaster) digunakan untuk mengidentifikasi kategori fungsional. Kategori fungsional didefinisikan dengan pengayaan setiap klaster oleh reaktom dan proses biologis GO. Untuk visualisasi gen resistensi di setiap jalur, semua gen di jalur diunggah ke mekanisme muatan STRING dan kode warna ditetapkan untuk lingkaran cahaya di sekitar simpul STRING.

2.11 Tes genetik untuk varian patogen
Genotipe untuk varian patogenik berulang (PV) pada BRCA1/2 dan gen kerentanan kanker tambahan dilakukan pada 17 pasien yang tidak memiliki informasi genetik sebelumnya. Panel gen menguji 51 varian, termasuk 30 PV BRCA1/2 dan 21 varian tambahan yang diketahui pada populasi Yahudi dan non-Yahudi. Pengujian dilakukan menggunakan sistem NanoCHIP XL (Savyon Diagnostics, Ashdod, Israel). Daftar lengkap varian yang diuji tersedia dalam referensi [ [ 47 ] ].

3 Hasil
3.1 Dataset longitudinal pasien kanker payudara yang menerima kemoterapi neoadjuvant
Sampel arsip yang cocok (tumor praperawatan, tumor pascaperawatan, dan epitel payudara normal yang berdekatan) dikumpulkan dari kohort pasien kanker payudara yang menjalani kemoterapi neoadjuvan (Gbr. 1A ). Yang penting, sampel diperkaya untuk sel tumor dengan makrodiseksi, dan sebagian besar sampel mencakup setidaknya 70% sel tumor, seperti yang dinilai oleh ahli patologi. Deskripsi lengkap karakteristik klinikopatologi dirangkum dalam Tabel 1 dan dirinci dalam Tabel S1 , termasuk selularitas tumor. Perawatan kemoterapi mencakup doksorubisin dan siklofosfamid berurutan standar perawatan diikuti oleh paclitaxel, dan Herceptin tambahan untuk pasien HER2+. Sampel pascaperawatan dinilai untuk respons oleh skor Miller–Payne (MP) [ [ 31 ] ] dan grade RCB [ [ 32 ] ]. Dataset mRNAseq seluruh transkriptom dibuat, menggunakan metodologi mRNAseq kami yang dioptimalkan, untuk sampel FFPE yang diarsipkan [ [ 30 ] ]. Diagram alir menyajikan kasus-kasus yang termasuk dalam analisis (Gbr. S1 ).

3.2 Berbagai jenis sampel dikelompokkan bersama yang menunjukkan status regulasi yang berbeda
Analisis komponen utama (PCA) berbasis ekspresi dari semua sampel menunjukkan pengelompokan sampel menurut jenis sampel (tumor pra- atau pasca-perawatan dan epitel normal) (Gbr. 1B ). Sampel pasca-perawatan diposisikan di antara jaringan normal yang berdekatan dan jaringan pra-perawatan kanker, yang mungkin mencerminkan rentang respons terhadap perawatan. Khususnya, payudara normal yang sehat dikelompokkan dengan sampel normal pasien kanker yang berdekatan (Gbr. 1B ), yang juga ditampilkan sebagai peta panas dari 1000 gen yang diekspresikan paling berbeda antara tumor dan normal (Gbr. S3 ). Pengamatan ini mendukung asumsi kami bahwa jaringan payudara normal yang berdekatan dari pasien kanker secara andal mewakili ekspresi jaringan payudara normal.

Modulasi keseluruhan yang terjadi selama perjalanan penyakit dan terapi diperkirakan dengan menghitung skor deregulasi jalur (PDS) untuk berbagai sampel, relatif terhadap semua sampel normal menggunakan Pathifier [ [ 42 , 43 ] ]. Algoritme ini menghitung PDS untuk setiap sampel, berdasarkan jarak dari proyeksi awan titik pada kurva utama (jalur representatif yang dicontohkan dalam Gambar 1C ). PDS dihitung untuk daftar gabungan jalur (proses GO_biological, KEGG, dan Reactome). Pandangan menyeluruh dari skor deregulasi ditunjukkan sebagai peta panas dari semua jalur di semua sampel, diurutkan berdasarkan pengelompokan hierarkis (Gambar 1D ). Khususnya, jaringan normal yang berdekatan dikelompokkan dengan jaringan payudara normal yang sehat, yang menunjukkan bahwa sampel normal pascaperawatan ini sebagian besar mewakili epitel payudara normal. Tumor praperawatan menunjukkan skor yang paling deregulasi, sedangkan tumor pascaperawatan bervariasi, beberapa dikelompokkan dengan sampel praperawatan, dan yang lainnya dengan normal yang berdekatan, terlepas dari selularitas tumor. Peta keseluruhan deregulasi jalur menunjukkan bahwa perubahan dinamis pada jalur tertentu bervariasi antara pasien. Untuk memeriksa divergensi dinamika jalur tertentu ini, kami memetakan PDS dari sampel yang cocok untuk setiap pasien (Gbr. 1E dan Gbr. S4 ).

Kami mengamati bahwa, sementara pada beberapa pasien, jalur kembali ke skor normal setelah perawatan, pada pasien lain jalur tetap pada keadaan deregulasi. Misalnya, jalur Anemia Fanconi sangat deregulasi pada pasien TNBC dan tetap deregulasi juga pascaperawatan, tetapi pada pasien ST-7 dan ST-1, PDS kembali ke tingkat normal pascaperawatan. Skor PDS memberikan pandangan global dari seluruh deregulasi jalur; namun, rute menuju resistensi adaptif dapat bervariasi di antara jalur dan pasien. Jadi, untuk memahami apakah variabilitas di antara pasien dalam jalur yang sama dikaitkan dengan gen tertentu atau gen yang beragam, kami melanjutkan ke analisis modulasi temporal pada tingkat gen per pasien.

3.3 Klasifikasi pola longitudinal untuk mengidentifikasi pola ekspresi gen yang terkait dengan respons
Dengan memeriksa modifikasi gen temporal, kami mengamati bahwa dinamika ekspresi gen, bukan tingkat absolutnya, sebagian besar bersifat spesifik pasien. Biasanya, gen menunjukkan dinamika ekspresi divergen yang berhubungan dengan respons (Gbr. 2A dan Gbr. S5 ) dengan pengecualian beberapa gen yang menunjukkan pola di seluruh pasien (Tabel S2 ; Gbr. S6A ). Beberapa gen dengan pola di seluruh pasien ini mengalami penurunan regulasi pada tumor primer relatif terhadap jaringan normal yang berdekatan, seperti yang divalidasi dalam kumpulan data TCGA, yang menunjukkan bahwa perbedaan ini bukan karena efek pengobatan pada jaringan normal yang berdekatan (Gbr. S6B ). Beberapa gen ini biasanya menunjukkan koekspresi pasangan gen tertinggi di antara pasien dalam kumpulan data titik waktu tunggal (TCGA dan METABRIC, Tabel S3 dan Gbr. S6C ), yang menunjukkan regulasi timbal baliknya.

Gambar 2
Buka di penampil gambar
Kekuatan Gambar
Klasifikasi pola longitudinal untuk mengidentifikasi pola ekspresi gen yang terkait dengan respons. (A) Sebagian besar gen menunjukkan dinamika ekspresi yang berbeda di antara pasien ( n  = 23) antara titik waktu, dengan pola yang biasanya terkait dengan respons pengobatan, seperti yang ditunjukkan dalam grafik garis. (B) Skema yang menyajikan hipotesis bahwa dinamika yang berbeda dalam ekspresi massal dapat menentukan gen yang terkait dengan resistensi terhadap pengobatan. (C) Gen yang diekspresikan secara berbeda diklasifikasikan ke dalam delapan kemungkinan pola modulasi teoritis. Berbagai pola dan skenario yang disarankan disajikan dalam warna yang berbeda. (D) Peta panas semua triplet gen per pasien setelah klasifikasi dengan mengkorelasikan nilai Transformasi Stabilisasi Varians dengan pola teoritis yang telah ditentukan sebelumnya. Setiap baris sesuai dengan satu gen pada pasien tertentu. Gen yang sama diplot untuk setiap pasien, menurut pola yang diklasifikasikan pada pasien individu. (E) Fraksi gen dalam berbagai pola untuk setiap pasien. Pasien diurutkan berdasarkan persentase gen dengan pola P7 dan P8, yang menunjukkan peningkatan disregulasi pascaperawatan relatif terhadap kadar normal. (F) Boxplot yang menunjukkan distribusi rasio “pola yang diatur ulang” terhadap “pola yang tidak diatur” (P1 + P2/P7 + P8), di tiga kelompok skor MP ( n  = 23 pasien) (garis hitam – rata-rata; kotak abu-abu muda – Standar; kotak abu-abu tua – 95% SEM). Rasio ini secara signifikan dikaitkan dengan skor respons MP (ANOVA satu arah, nilai- P = 0,0095). (G) Kurva Kaplan–Meier dari kelangsungan hidup bebas kekambuhan untuk seluruh kohort ( n  = 23 pasien) berdasarkan rasio fraksinya: (P1 + P2/P7 + P8), titik potong = 5, nilai- P log rank 0,09. MP, skor respons Miller–Payne.
Untuk sebagian besar gen, pola temporal, bukan nilai ekspresi absolut, berkorelasi dengan skor respons pengobatan. Misalnya, tingkat ekspresi praperawatan absolut FANCA atau KIT (Gbr. 2A dan Gbr. S5 ) tidak terpisah antara berbagai kelompok respons MP patologis. Sebaliknya, pola ekspresi dinamis dibedakan antara kelompok respons (Gbr. 2A dan Gbr. S5 ). Divergensi dalam dinamika ini membuat kami berhipotesis bahwa mengklasifikasikan gen berdasarkan pola ekspresinya dapat mengidentifikasi gen terkait resistensi, seperti yang digambarkan dalam skema pada Gbr. 2B . Dengan demikian, kami melakukan analisis pola dinamis untuk gen yang diekspresikan secara berbeda, mirip dengan penelitian kami sebelumnya [ [ 28 , 29 ] ], memanfaatkan kekuatan kumpulan data longitudinal. Kami mendefinisikan delapan kemungkinan pola teoritis (P1–P8) dari modulasi ekspresi gen melalui triplet yang cocok dari tahap pengobatan tumor normal (Gbr. 2C ). Setiap pola mewakili skenario kejadian yang berbeda melalui perkembangan tumor dan tahap pengobatan:
Tingkat ekspresi yang meningkat (P1) atau menurun (P2) dari sampel normal ke sampel yang telah diobati sebelumnya dan kembali ke tingkat normal setelah pengobatan. Hal ini mungkin merupakan reregulasi gen setelah pengobatan.
Tingkat ekspresi yang meningkat (P3) atau menurun (P4) dari normal pada sampel yang telah diberi perlakuan awal, yang tetap serupa dengan perlakuan akhir, mungkin menunjukkan keadaan resistensi.
Ekspresi serupa antara sampel normal dan tumor, yang menurun (P5), atau meningkat (P6) hanya setelah pascaperawatan, yang dapat disimpulkan sebagai efek perawatan.
Peningkatan konstan (P7) atau penurunan konstan (P8) dalam tingkat ekspresi dari normal melalui sampel yang telah diobati sebelumnya ke sampel setelah diobati. Hal ini dapat menunjukkan memburuknya proses kanker.
Lintasan gen per pasien diklasifikasikan ke dalam berbagai pola dengan mengorelasikan nilai Variance Stabilizing Transformation (VST) (contoh pola teoritis dan gen yang ditetapkan pada Gambar S7 ). Peta panas yang menyajikan distribusi di seluruh pola semua lintasan gen per pasien ditunjukkan pada Gambar 2D . Gen yang sama diplot untuk semua pasien, menurut pola klasifikasi spesifik pasien. Untuk memvalidasi arah pola klasifikasi, kami membandingkan perubahan lipatan antara pasangan tumor-normal dengan perubahan lipatan yang diamati dalam kumpulan data METABRIC (perubahan lipatan 1,5 dan tingkat penemuan palsu [FDR] ≤ 0,05) [ [ 8 ] ] dan menemukan kesepakatan dalam arah peningkatan regulasi atau penurunan regulasi pada tumor vs jaringan normal.

3.4 Fraksi pola gen longitudinal berhubungan dengan respons
Untuk memeriksa hubungan antara pola ekspresi dan respons, kami memetakan fraksi gen yang ditetapkan untuk setiap pola per pasien (Gbr. 2E ). Pasien dengan skor respons yang lebih baik (nilai MP yang lebih tinggi) menunjukkan fraksi gen yang lebih tinggi dengan pola P1 dan P2, di mana ekspresi kembali ke tingkat normal. Responden yang buruk menunjukkan fraksi P3 atau P4 yang lebih tinggi dan terutama P7 atau P8, yang menunjukkan peningkatan dalam proses terkait kanker pascaperawatan. Khususnya, selularitas sel kanker dalam sampel tidak terkait dengan fraksi pola. Kami menemukan bahwa rasio antara fraksi gen dengan pola P1 dan P2 relatif terhadap P7 dan P8 (P1 + P2/P7 + P8), secara signifikan terkait dengan skor respons MP (analisis varians [ANOVA] uji, P -nilai = 0,0095) (Gbr. 2F ). Fraksi pola ini juga terkait dengan kelangsungan hidup bebas kekambuhan, meskipun tidak signifikan karena ukuran kohort yang kecil (Gbr. 2G ). Skor MP untuk kelompok ini juga dikaitkan dengan kelangsungan hidup bebas kekambuhan (Gbr. S8 , ns), yang menekankan validitas klinisnya [ [ 31 ] ]. Secara keseluruhan, data menunjukkan bahwa fraksi gen yang lebih tinggi dengan pola yang diatur ulang dikaitkan dengan respons yang lebih baik dan fraksi gen yang lebih tinggi dengan pola yang tidak diatur dikaitkan dengan respons yang buruk. Temuan ini memotivasi analisis pola ekspresi untuk mengidentifikasi gen yang terkait dengan resistensi.

3.5 Analisis pola gen mengidentifikasi gen dan jalur yang terkait dengan respons yang menyatu menjadi modul resistensi terbatas
Selanjutnya kami berusaha untuk menentukan gen yang mungkin memiliki peran dalam resistensi adaptif dengan mencari gen dengan pola disregulasi berulang di seluruh pasien. Untuk memahami bagaimana berbagai pola didistribusikan di antara pasien, kami mengelompokkan matriks jumlah pasien per pola untuk semua gen yang diekspresikan secara berbeda (~6900 gen). Kami mengamati bahwa, seringkali, gen diklasifikasikan menjadi satu pola utama atau pola bimodal di seluruh pasien (Gbr. 3A ). Kami berhipotesis bahwa suatu gen dapat dikaitkan dengan resistensi adaptif jika pola utamanya di sebagian besar pasien (>50%) adalah salah satu pola disregulasi P3/P7 atau P4/P8, yaitu, bahwa ekspresi gen terkait kanker bertahan atau meningkat setelah pengobatan. Kami menemukan 198 gen dengan pola disregulasi persisten yang berulang di seluruh pasien (Tabel S4 ). Sebagai alternatif, gen yang terkait dengan resistensi juga dapat menunjukkan pola bimodal di seluruh pasien, di mana gen yang sama dapat disregulasi pada satu pasien dan diregulasi ulang pada pasien lain. Kami berspekulasi bahwa distribusi bimodal di antara pasien ini dapat membedakan antara responden yang baik dan responden yang buruk (Gbr. 3B ). Dengan demikian, kami membandingkan antara skor respons MP dari responden yang baik dan responden yang buruk untuk gen dengan pola bimodal menggunakan uji jumlah Wilcoxon Rank dan mengidentifikasi 254 gen signifikan (Gbr. 3B ) (Tabel S5 ). Untuk memperhitungkan bias subtipe, kami mengulangi uji ini hanya untuk pasien HR+ (20 pasien) dan menemukan 122 gen signifikan (Tabel S5 ).

Gambar 3
Buka di penampil gambar
Kekuatan Gambar
Analisis pola gen mengidentifikasi gen dan jalur yang terkait dengan resistensi yang konvergen menjadi modul resistensi terbatas. (A) Bagan alir algoritma analisis pola gen: Gen yang diekspresikan secara diferensial yang diklasifikasikan ke pola dinamis biasanya menunjukkan satu pola utama atau pola bimodal, (misalnya, P1 atau P3). (B) Untuk gen dengan pola bimodal, kami membandingkan skor respons patologis MP pasien (menggunakan uji jumlah Wilcoxon Rank). Kami mengidentifikasi 254 gen yang polanya berbeda secara signifikan antara responden yang baik dan responden yang buruk. (C) Analisis jaringan (STRING) untuk semua 452 gen yang terkait dengan resistensi (198 pola utama yang bertahan dan 254 gen yang membedakan responden yang baik dan yang buruk) menunjukkan konvergensi gen resistensi menjadi sejumlah modul resistensi yang terbatas. Yang penting, tiga modul sangat terkait dengan mekanisme kerja kemoterapi yang diberikan.
Untuk menilai apakah pengobatan mempengaruhi sampel normal yang berdekatan, kami membandingkan tingkat ekspresi 452 gen terkait resistensi di seluruh sampel normal (baik yang berdekatan maupun yang sehat) dan tumor pascapengobatan. Analisis kami mengungkapkan bahwa jaringan payudara normal yang berdekatan dari kelompok pasien kanker dengan jaringan payudara normal yang sehat secara jelas terpisah dari tumor pascapengobatan. Pola pengelompokan ini menunjukkan efek pengobatan minimal pada jaringan normal yang berdekatan (Gbr. S9A ). Kami lebih jauh menekankan temuan ini dengan berfokus pada pasien dengan kanker payudara bilateral (ST-52), di mana pemisahan antara sampel tumor normal yang berdekatan dan pascapengobatan sangat jelas (Gbr. S9B ). Untuk memahami konektivitas fungsional antara gen terkait resistensi yang teridentifikasi (total 452), kami melakukan analisis jaringan STRING. Kami menemukan bahwa gen-gen tersebut konvergen ke sejumlah modul resistensi yang terbatas (Gbr. 3C ). Khususnya, tiga modul terkait dengan mekanisme kerja kemoterapi yang diberikan (taxane, cyclophosphamide dan doxorubicin) [ [ 48 ] ]. Taxanes bekerja dengan mengganggu fungsi mikrotubulus, sehingga menghambat mitosis; doksorubisin menargetkan replikasi DNA dan siklofosfamid memberikan efek sitotoksiknya terutama dengan mengikat silang untaian DNA. Kami mengidentifikasi tujuh gen yang diperkaya dalam jalur perbaikan DNA, terutama dalam jalur Anemia Fanconi (Tabel S6 ). Modul resistensi lain yang muncul terkait dengan remodeling matriks ekstraseluler (ECM), pensinyalan PI3K-AKT dan MAPK, regulasi angiogenesis, transpor glukosa, dan metabolisme lipid (Gbr. 3C ).

Bahasa Indonesia: Untuk lebih jauh mengeksplorasi konteks fungsional gen-gen yang berhubungan dengan resistensi, kami melakukan analisis pengayaan set gen dengan clusterprofiler [ [ 45 ] ] menggunakan set data gabungan (proses GO_biological KEGG dan Reactome) yang diunduh dari EnrichR [ [ 44 ] ]. Kami menghitung pengayaan secara terpisah untuk setiap set pola bimodal (P1_P3_P7 dan P2_P4_P8) dan mengidentifikasi 152 jalur yang secara signifikan berhubungan dengan resistensi ( nilai q < 0,05) (lihat Tabel S6 ). Kami kemudian mengulangi pengayaan untuk gen-gen yang berhubungan dengan resistensi yang signifikan untuk pasien HR+, mengidentifikasi 78 jalur yang diperkaya ( nilai q < 0,05) (Tabel S6 ). Untuk mendapatkan daftar jalur yang diperkaya yang kuat, dengan memperhitungkan heterogenitas subtipe, kami menggunakan interseksi antara kedua analisis. Jalur tersebut terkait dengan siklus sel, mitosis, replikasi DNA dan perbaikan DNA serta kaskade MAPK dan PI3K, transpor glukosa, transpor lipid, dan organisasi matriks ekstraseluler.

3.6 Setiap pasien menunjukkan disregulasi yang berbeda dari berbagai gen dalam jalur
Untuk memahami divergensi gen yang mengalami disregulasi antara pasien, kami berfokus pada jalur yang membedakan responder baik dan buruk dan memeriksa semua gen yang diekspresikan dalam jalur yang mengalami disregulasi, untuk setiap pasien individu. Kami mendefinisikan gen dengan pola P1/P2 sebagai “gen yang diatur ulang” dan gen dengan pola P3/P4/P7/P8 sebagai “gen resistensi.” Kami mengamati bahwa fraksi dan identitas gen resistensi dalam jalur yang sama bervariasi antara pasien. Oleh karena itu, kami menghitung skor resistensi untuk setiap pasien—fraksi gen resistensi per jalur. Signifikansi skor ini dihitung dengan melakukan bootstrapping pada set gen acak dalam ukuran jalur yang sama pada setiap pasien (diambil sampelnya 10.000 kali). Perhitungan skor resistensi untuk dua pasien dicontohkan oleh peta jaringan jalur representatif (Reactome: Replikasi DNA), dan nilai- P diwakili oleh histogram fraksi acak gen resistensi untuk jalur ini (Gbr. 4A ).

Gambar 4
Buka di penampil gambar
Kekuatan Gambar
Skor resistensi jalur dan peta kemoresistom individual menunjukkan divergensi antar pasien. (A) Untuk setiap jalur terkait resistensi pada setiap pasien, gen diberi mode resistensi atau reregulasi berdasarkan pola klasifikasinya. Kami kemudian menghitung skor resistensi: fraksi “gen resistan” dalam jalur dan nilai – P untuk skor ini, dihitung dengan melakukan bootstrapping pada sekumpulan gen acak dalam ukuran jalur yang sama (diambil sampelnya 10.000 kali). Contoh skor resistensi terhitung pada dua pasien representatif untuk jalur replikasi DNA, disajikan sebagai jaringan. Titik merah: gen resistensi, titik hijau: gen yang diregulasi ulang, titik abu-abu: gen tanpa pola yang ditetapkan. Perhitungan nilai- P disajikan sebagai histogram distribusi fraksi acak pada pasien ini. Garis merah menunjukkan skor resistensi (fraksi gen resistan dalam jalur ini untuk pasien ini). (B) Skor resistensi dan nilai- P untuk 24 jalur representatif yang diperkaya untuk gen terkait resistensi disajikan sebagai diagram balon. Setiap baris mewakili jalur dan setiap kolom mewakili pasien. Ukuran lingkaran menunjukkan skor resistensi dan kode warna menunjukkan nilai – P dari skor resistensi. Parameter klinis untuk setiap pasien diberi kode warna di panel atas. Jalur dikumpulkan menurut modul fungsionalnya, dengan warna yang sama seperti pada jaringan STRING di A. (C) Peta kemoresistom personal untuk pasien representatif yang menggambarkan divergensi resistensi adaptif untuk setiap pasien. Peta tersebut mewakili kategori resistensi utama, dihitung untuk jalur yang paling signifikan di setiap kategori. Skor resistensi dilambangkan dengan tinggi batang dan nilai- P diberi kode warna. HP, hormon positif; RCB, beban kanker residual; TN, triple-negatif; TP, triple-positif.
Skor resistensi dan nilai- P per pasien untuk daftar 24 jalur terkait resistensi yang signifikan yang dikurasi disajikan oleh plot balon, yang menunjukkan variabilitas dalam jalur yang tidak teratur antara pasien (Gbr. 4B dan Tabel S7 ). Yang penting, terbukti dari plot bahwa pada beberapa pasien beberapa jalur menyebabkan resistensi, sementara pada pasien lain hanya beberapa jalur yang menunjukkan skor resistensi yang tinggi. Perbedaan yang jelas antara subtipe terlihat pada pasien TNBC (ST_42, ST_57, dan ST_60), yang menunjukkan beberapa jalur terkait resistensi yang biasanya terkait dengan pembelahan sel, sementara pada pasien HR+, jalur terkait resistensi bervariasi dan spesifik untuk setiap pasien. Meskipun variasi antara subtipe nyata, variabilitas antar pasien juga terlihat antara pasien yang memiliki subtipe yang sama dan menerima perawatan yang sama. Misalnya, pasien ST-11 dan pasien ST-35 menunjukkan skor resistensi yang signifikan hanya dalam beberapa jalur, yaitu jalur Anemia Fanconi pada ST-11 dan jalur transporter ABC pada ST-35 (Gbr. 4B ). Sebaliknya, pasien ST-10, ST-32, dan ST-7 menunjukkan skor resistensi yang signifikan dalam beberapa jalur. Meskipun hanya tiga pasien yang memiliki TNBC, penyertaan subtipe ini dalam kelompok tersebut memungkinkan penentuan jalur yang kurang terwakili dalam HR+ tetapi terbukti pada beberapa pasien HR+, seperti pembelahan mitosis pada ST-1, organisasi mikrotubulus pada ST-32, dan perbaikan DNA pada ST-8 (Gbr. 4B ). Divergensi dalam jalur disregulasi spesifik ini sebagai respons terhadap kemoterapi yang sama menunjukkan kemunculan resistensi adaptif secara individu.

Untuk menyelidiki hubungan potensial antara resistensi yang muncul dalam jalur perbaikan DNA dan mutasi spesifik yang terkait dengan defisiensi rekombinasi homolog, kami melakukan pengujian genetik menggunakan panel 50 mutasi. Panel ini mencakup 30 mutasi BRCA1/2, yang dipilih untuk mencakup varian yang paling umum dalam populasi studi kami berdasarkan etnis. Di antara 22/23 pasien yang diuji, 20 negatif untuk semua mutasi yang disaring. Tiga pasien menunjukkan hasil positif: ST-1 dan ST-28 mengandung mutasi CHEK2, sementara ST-52 membawa mutasi BRCA2 (Tabel S1 ). Khususnya, skor resistensi yang signifikan untuk jalur Anemia Fanconi atau jalur perbaikan ikatan silang interstrand diamati pada pasien ST-11, ST-42, ST-57, dan ST-28. Analisis kami tidak menemukan hubungan yang jelas antara mutasi genetik yang teridentifikasi dan skor resistensi yang meningkat dalam jalur ini. Temuan ini menunjukkan bahwa resistensi adaptif yang diamati dalam jalur perbaikan DNA mungkin tidak secara langsung disebabkan oleh mutasi yang sudah ada sebelumnya pada gen yang kami uji.

Untuk menyorot kontribusi berbeda dari setiap modul terhadap munculnya resistensi pada masing-masing pasien, kami membuat peta kemoresistom pribadi dengan diagram batang melingkar (Gbr. 4C ). Peta tersebut menggambarkan koeksistensi beberapa modul/mekanisme resistensi pada masing-masing pasien, sambil menyorot variabilitas antar pasien melalui dampak masing-masing modul. Yang penting, beberapa modul menunjukkan kejadian spesifik pasien atau subtipe, yang menekankan munculnya resistensi secara individu. Kami mengamati bahwa modul resistensi yang terkait dengan siklus sel dan pembelahan sel berlimpah pada tumor TNBC dan tumor yang membelah cepat, sementara modul lain, seperti remodeling ECM dan regulasi PI3K/AKT, lebih merupakan karakteristik tumor HR+. Meskipun demikian, modul proliferasi sel juga terlihat pada pasien non-TNBC, tetapi seringkali lebih sporadis dan spesifik terhadap mekanisme kerja satu obat.

3.7 Variabilitas dalam disregulasi gen dan jalur mengungkapkan prinsip-prinsip dalam resistensi kemoterapi adaptif
Untuk lebih jauh mengeksplorasi variabilitas antar pasien dalam jalur yang tidak teratur, kami meneliti secara mendalam berbagai gen dalam setiap jalur yang terkait dengan resistensi. Ketika membandingkan gen resistensi antar pasien, kami menemukan bahwa terdapat pusat gen resistensi (Pola P3, P4, P7, dan P8), yang muncul kembali pada banyak pasien. Pusat-pusat ini direpresentasikan sebagai peta panas untuk jumlah pasien dengan pola resistensi per gen (Gbr. 5A: panel kanan dan Gbr. S10 ).

Gambar 5
Buka di penampil gambar
Kekuatan Gambar
Variabilitas dalam penataan ulang gen dan disregulasi jalur mengungkap prinsip resistensi adaptif setelah kemoterapi. (A) Contoh variabilitas dalam disregulasi jalur antara pasien. Presentasi jaringan gen yang diekspresikan dalam dua jalur representatif, diwarnai oleh analisis pola gen untuk mengklasifikasikan gen menjadi gen yang diregulasi ulang atau gen resistensi (lihat Gambar 4B ). Sementara pasien ST_42 (TNBC) menunjukkan beberapa gen resistensi di kedua jalur, pasien ST_7 dan pasien ST_52 menunjukkan kecenderungan resistensi yang berlawanan untuk kedua jalur tersebut. Di sebelah kanan: jalur yang sama, direpresentasikan sebagai peta panas untuk jumlah pasien dengan pola resistensi per gen, yang menunjukkan hub gen untuk resistensi di jalur tersebut. (B) Prinsip resistensi adaptif setelah kemoterapi. Ratusan gen mengalami penataan ulang setelah kemoterapi, dengan cara yang spesifik untuk pasien. Gen yang terkait dengan resistensi sering kali disiapkan secara transkripsi dalam tumor primer. Jalur yang sama dapat mengalami disregulasi dengan menata ulang gen yang berbeda. Gen yang diubah dan jalur yang tidak terkendali berkumpul menjadi modul resistensi yang terbatas, biasanya terkait dengan mekanisme kerja obat yang diberikan.
Variabilitas antar pasien dari gen yang diubah susunannya dicontohkan untuk dua jalur oleh peta jaringan gen (Gbr. 5A ). Misalnya, dalam jalur Anemia Fanconi, banyak gen menunjukkan pola resistensi pada pasien ST_42 (TNBC), sementara pada pasien ST_7 (HR+) dan ST_52 (HR+) jumlah gen resistensi bervariasi (Gbr. 5A ). Secara khusus, gen dalam jalur ini untuk pasien ST_7 menunjukkan sebagian besar pola yang diatur ulang, meskipun semua pasien menerima siklofosfamid. Intervariabilitas antara pasien ini juga terbukti dalam jalur depolymerisasi mikrotubulus. Misalnya, dalam jalur Anemia Fanconi, RMI2 dan FANCA menunjukkan pola resistensi pada 13 pasien (Gbr. S11A,B ). Yang penting, pasien ST_52 tidak menerima paclitaxel dan kami hanya mengidentifikasi gen yang diatur ulang dalam jalur depolymerisasi mikrotubulus (Gbr. S11C ), yang menunjukkan bahwa pengaturan ulang terkait resistensi diinduksi oleh Taxane, yang menargetkan jalur ini. Namun, variabilitas antar pasien tinggi dan kami hanya mengamati gen yang diatur ulang pada pasien lain (kebanyakan HR+/HER2+) yang menerima paclitaxel (Gbr. S11C ).

Secara keseluruhan, hasil penelitian ini mengarahkan kami untuk menyarankan prinsip resistensi adaptif pascakemoterapi (diilustrasikan dalam Gambar 5B ). Kami mengusulkan bahwa ratusan gen mengalami penataan ulang setelah kemoterapi. Gen yang ditata ulang sebagian besar adalah gen terkait kanker yang telah dimodulasi selama proses karsinogenik relatif terhadap jaringan epitel normal, beberapa di antaranya kembali ke tingkat normal, sementara yang lain mempertahankan tingkat yang tidak teratur. Meskipun banyak gen dan jalur menunjukkan disregulasi pascaperawatan, mereka bertemu pada jumlah modul resistensi terbatas yang terbatas. Lebih jauh lagi, gen yang ditata ulang dan jalur yang tidak teratur bersifat spesifik pasien; dengan demikian, gen yang berbeda dalam jalur yang sama dapat ditata ulang pada pasien individu.

4 Diskusi
Hasil penelitian ini memberikan perspektif yang luas tentang kompleksitas adaptasi sel kanker terhadap kemoterapi. Desain penelitian dan pendekatan analisis kami berbeda dari penelitian sebelumnya dalam dua hal: (a) Memeriksa perubahan yang berhubungan dengan epitel normal yang sesuai; (b) Mendekonstruksi jalur menuju modulasi pola dinamis untuk setiap gen per pasien, dengan demikian mengungkap berbagai rute munculnya resistensi adaptif pada masing-masing pasien. Memeriksa modulasi spesifik pasien sebelum dan sesudah pengobatan relatif terhadap epitel normal dari pasien yang sama mengungkapkan bahwa pengkabelan ulang ekspresi gen adalah interaksi komposit, unik untuk setiap tumor. Kami mendalilkan bahwa berbagai jalur, yang sudah mengalami disregulasi pada tumor primer, mempertahankan keadaan toleran obat dengan mengaktifkan banyak gen yang menjaga aktivitas abnormal mereka. Dengan demikian, untuk bertahan hidup dari aktivitas toksik suatu obat, sel tumor mempertahankan keadaan disregulasi atau mengintensifkannya, khususnya melewati gangguan obat.

Pendekatan kami berupa pengambilan sampel longitudinal, dikombinasikan dengan analisis pola, menekankan pentingnya memeriksa perubahan dinamis untuk mengungkap mekanisme resistensi obat. Sebagian besar gen resistensi tidak dapat dideteksi dalam analisis titik waktu tunggal, karena kadar absolutnya pada tumor primer tidak dapat dibedakan antara responden yang baik dan yang buruk. Misalnya, gen TOP2A, yang secara spesifik ditargetkan oleh doksorubisin, sebelumnya dievaluasi sebagai penanda manfaat dari obat ini [ [ 49 ] ], tetapi hasil yang tidak konsisten mengganggu kemanjuran prediktifnya. Data kami dengan jelas menunjukkan bahwa kadar ekspresi TOP2A pada tumor primer tidak konsisten dengan respons, tetapi pola ekspresi longitudinalnya membedakan antara responden yang baik dan yang buruk (Gbr. S5 ). Ini mungkin menjelaskan mengapa fenotipe tumor heterogen menantang identifikasi penanda molekuler dan menunjukkan bahwa perubahan dinamis spesifik pasien sangat penting untuk menargetkan gen dan jalur resistensi.

Keuntungan utama dari kumpulan data kami adalah penyertaan jaringan normal berdekatan yang cocok. Dalam studi longitudinal sebelumnya, pencarian gen yang diekspresikan secara berbeda sebelum dan sesudah pengobatan melewatkan gen yang tidak dimodifikasi yang mempertahankan keadaan disregulasi, seperti pada Pola P3. Di sisi lain, gen yang sering dimodifikasi pascapengobatan tetapi sebenarnya diregulasi ulang (pola P1/P2), dapat disalahartikan sebagai penanda resistensi kemoterapi (dicontohkan dalam skema, Gambar S12 ). Misalnya, CYR61 sebelumnya dilaporkan mengalami peningkatan regulasi pasca-NAT [ [ 23 ] ], dan disarankan sebagai penanda resistensi kemoterapi. Dalam kumpulan data kami, kadarnya memang meningkat regulasi pasca-NAT, tetapi pada kadar yang sama seperti pada jaringan normal berdekatan, sehingga disimpulkan sebagai gen yang diregulasi ulang (Gambar S12A ). Sebaliknya, upregulasi FN1 pasca-NAT sebelumnya dikaitkan dengan respons yang buruk [ [ 50 ] ] dan terutama menunjukkan pola resistensi dalam kumpulan data kami (Gbr. S12B ). Yang penting, jaringan normal yang berdekatan dapat dipengaruhi oleh proses karsinogenik di sekitarnya, sebuah fenomena yang dikenal sebagai kankerisasi lapangan, atau mungkin cenderung terkena kanker karena faktor genetik. Oleh karena itu, sangat penting untuk memvalidasi bahwa jaringan normal yang berdekatan mirip dengan jaringan payudara normal yang sehat, seperti yang ditunjukkan pada Gbr. S3 . Namun, dalam studi longitudinal yang membandingkan setiap tumor dengan normal yang berdekatan adalah referensi yang paling relevan. Demikian pula, penting untuk memastikan bahwa perubahan yang diamati pada tumor pascaperawatan khusus untuk jaringan kanker dan tidak secara global mempengaruhi jaringan normal yang berdekatan, seperti yang diilustrasikan pada Gbr. S9 .

Fakta bahwa sebagian besar gen resistensi diidentifikasi dalam pola P3/P7 atau P4/P8, dimodulasi dalam keadaan karsinogenik relatif terhadap epitel normal, menunjukkan bahwa tumor resisten terhadap pengobatan dengan mempertahankan jalur yang tidak teratur. Jalur pensinyalan yang penting untuk kelangsungan hidup sel kanker mempertahankan aktivitasnya atau beralih ke gigi yang lebih tinggi untuk bertahan hidup dalam kondisi kemoterapi yang ekstrem. Pengamatan kami menunjukkan bahwa subpopulasi persisten memiliki profil transkripsi yang disiapkan, relatif terhadap keadaan normal, yang menyampaikan bias ekspresi gen terhadap kelangsungan hidup. Gagasan ini sesuai dengan studi in vitro terbaru dari sel-sel yang toleran terhadap obat, yang menunjukkan bahwa resistensi fenotipik ditentukan apriori dengan cara nonstokastik dalam subpopulasi sel persister [ [ 51 , 52 ] ]. Studi lain tentang sekuensing sel tunggal pasien TNBC sebagai respons terhadap kemoterapi mendukung konsep ini, yang menunjukkan seleksi klonal dari genotipe yang sudah ada sebelumnya, serta pemrograman ulang transkripsi tanda-tanda resisten dalam klon persisten [ [ 6 ] ]. Mirip dengan hasil kami, mereka menemukan bahwa program transkripsi menyatu ke beberapa jalur resistensi dan sebagian besar diperoleh setelah pengobatan. Data kami mengonfirmasi konsep bahwa fenotipe resistan sudah diperoleh pada tumor primer, relatif terhadap jaringan normal yang cocok, dan disregulasinya telah bertahan atau meningkat. Namun, penting untuk menekankan bahwa meskipun gen-gen ini disiapkan pada tumor primer, pada banyak pasien ekspresinya dinormalisasi, sehingga membuat profil tumor primer tidak cukup untuk menentukan gen-gen ini atau jalur hilirnya sebagai mekanisme resistensi yang potensial.

Perubahan dinamis dalam respons terhadap kemoterapi dapat terjadi sebagai hasil dari seleksi genetik atau oleh adaptasi transkripsi besar-besaran. Kinetika adaptasi transkripsi jauh lebih cepat daripada evolusi genetik, dan diasumsikan bahwa mekanisme nongenetik sementara pada akhirnya akan tertanam melalui mekanisme epigenetik atau genetik jangka panjang [ [ 2 ] ]. Analisis pola in vitro dari omik waktu-waktu selama pengobatan mengungkapkan adaptasi transkripsi besar-besaran segera diikuti oleh perubahan epigenetik yang sesuai, disarankan untuk menstabilkan fenotipe yang resistan [ [ 53 ] ]. Khususnya, penelitian kami mengukur modulasi dalam tingkat ekspresi massal yang dapat mewakili seleksi transkripsi sel yang mengekspresikan gen atau adaptasi transkripsi dengan penataan ulang gen. Apa pun itu, tingkat bersih menunjukkan keadaan disregulasi atau penataan ulang tumor. Seperti yang ditunjukkan dalam penelitian kami dan penelitian sebelumnya [ [ 54 ] ], penataan ulang gen menyatu menjadi serangkaian fenotipe yang resistan terhadap obat yang terbatas. Akan tetapi, sementara sebagian besar penelitian berfokus pada gen dan jalur yang dimiliki bersama di antara pasien, kami menekankan heterogenitas fenotipik antara pasien dan perbedaan rute resistensi adaptif.

Munculnya jalur pensinyalan kunci dalam kumpulan data kami, yang secara langsung terkait dengan mekanisme kerja kemoterapi, menyerupai mekanisme resistensi terhadap inhibitor kinase, yang dengannya aktivasi target obat yang terus-menerus melewati efek penghambatannya [ [ 1 ] ]. Meskipun demikian, hasilnya juga menunjukkan perolehan mekanisme resistensi tambahan, seperti angiogenesis, terkait matriks ekstraseluler (ECM), pensinyalan PI3K/AKT, dan perubahan metabolisme lipid. Sementara pada tumor TNBC jalur proliferasi sebagian besar meningkat (P1), pada tumor HR+ jalur diferensiasi menurun (P2). Pada tumor TNBC jalur resistensi dominan terkait dengan siklus sel, spindel mitosis, jalur Hedgehog, replikasi DNA, dan perbaikan DNA, yang menunjukkan tumor yang sangat proliferatif ini. Pada tumor HR+, kami terutama mengamati pola resistensi pada jalur PI3K/Akt, regulasi kaskade MAPK, jalur terkait ECM, transpor lipid, dan jalur transporter ABC. Memang, aktivitas nongenomik estrogen telah disarankan untuk menginduksi resistensi kemoterapi, mengaktifkan PI3K/AKT, dan memodifikasi respons kerusakan DNA [ [ 55 ] ]. Pola resistensi muncul pada subkelompok spesifik transporter ABC, keluarga ABCA2, yang terlibat dalam transportasi steroid yang dapat mengatur proliferasi dan kelangsungan hidup kanker payudara, seperti yang disarankan dalam studi kanker prostat [ [ 56 ] ].

Temuan kami mungkin memiliki beberapa implikasi klinis di masa mendatang. Mengatasi kondisi toleran obat dalam pengaturan klinis merupakan tantangan, mengingat tidak ada satu mutasi pendorong untuk ditargetkan. Data kami menunjukkan bahwa sementara banyak fenotipe resistan muncul secara paralel, jalur yang diatur ulang mengungkapkan sensitivitasnya, memungkinkan kami untuk secara rasional memilih rencana terapi berikutnya. Tumor residual setelah NAT memiliki tingkat kejadian yang lebih tinggi untuk kekambuhan [ [ 15 ] ]. Tumor residual ini memiliki profil molekuler yang berbeda yang mengungkapkan resistensi yang diperolehnya di satu sisi, dan kerentanannya (yaitu, gen yang diatur ulang setelah pengobatan) di sisi lain. Oleh karena itu, pembuatan profil kemoresistom segera setelah operasi dapat memandu rencana pengobatan adjuvan dari terapi yang ditargetkan untuk mencegah penyakit berulang atau metastasis. Meskipun pembuatan profil biopsi pra dan pasca pengobatan dalam praktik klinis merupakan tantangan, perkembangan masa depan seperti RNA bebas sel atau pembuatan profil epigenetik DNA tumor dapat mengarah pada pembuatan profil kemoresistom selama siklus awal pengobatan neoadjuvan. Selain itu, kemajuan terkini memungkinkan untuk mengevaluasi berbagai perawatan ex vivo pada biopsi tumor yang berdekatan dengan diagnosis [ [ 57 ] ]. Pemetaan kemoresistom berdasarkan pola ekspresi gen sampel tumor yang dirawat ex vivo sebelum NAT dapat berfungsi untuk mengoptimalkan rencana perawatan yang dipersonalisasi. Perawatan ex vivo memungkinkan estimasi resistensi terhadap agen tunggal dan kombinasi obat. Untuk efektivitas biaya, membuat profil gen kemoresistom daripada seluruh transkriptom mungkin bermanfaat, meskipun dapat membatasi luasnya informasi. Misalnya, mengidentifikasi skor resistensi yang tinggi dalam jalur depolymerisasi mikrotubulus dapat mengindikasikan potensi resistensi adaptif untuk taxanes, sementara skor resistensi yang tinggi pada gen Anemia Fanconi dapat mengindikasikan munculnya resistensi adaptif terhadap siklofosfamid. Pendekatan ini dapat diperluas untuk membuat peta resistom untuk imunoterapi (imunoresistom) atau terapi bertarget lainnya, yang memungkinkan pemilihan perawatan yang lebih tepat.

Hasil studi kami dibatasi oleh ukuran kohort yang kecil dan subtipe yang berbeda-beda. Temuan dan prinsip yang disarankan harus divalidasi lebih lanjut dalam kohort yang lebih besar, serta dalam rejimen pengobatan lainnya. Sementara ukuran kohort dapat membatasi generalisasi, kerangka kerja longitudinal memberikan wawasan berharga tentang mekanisme resistensi khusus pasien, yang merupakan langkah penting menuju pemahaman dan kemungkinan mengatasi resistensi adaptif dalam konteks pengobatan yang dipersonalisasi. Meskipun ukuran kohort terlalu kecil untuk analisis subkelompok, kami mengklaim bahwa analisis subtipe gabungan menguntungkan untuk mendeteksi modulasi khusus pasien. Banyak modulasi terkait resistensi yang spesifik subtipe, terutama di TNBC, ditemukan pada beberapa pasien dari subtipe lain. Selain itu, studi ini difokuskan pada NAT, di mana pengobatan sistemik diberikan dengan adanya tumor primer. Sulit untuk menentukan apakah konsep resistensi adaptif, seperti yang diperoleh dari kumpulan data NAT, berlaku untuk pengaturan pengobatan adjuvan, di mana tidak ada tumor yang terbukti secara klinis. Terakhir, meskipun kami memperkaya sampel kami untuk sel epitel dengan pembedahan makro, faktor pengganggu yang dapat mengganggu kesimpulan analisis kami terkait dengan perubahan terkait stroma setelah perawatan, termasuk representasi berlebihan sel kolagen dan sel imun. Keterbatasan signifikan data RNA-seq tumor massal adalah bahwa data tersebut mewakili campuran berbagai jenis sel. Sementara data ekspresi gen kami dikelompokkan menurut jenis sampel dan tidak secara terang-terangan dipengaruhi oleh selularitas tumor, kami tidak dapat sepenuhnya mengesampingkan kemungkinan bahwa beberapa pola ekspresi gen yang diamati mencerminkan modulasi pada jenis sel lain atau perubahan proporsinya pascaperawatan. Selain itu, keterbatasan lainnya adalah potensi representasi berlebihan subklon dominan dari domain spasial yang diambil sampelnya, khususnya subklon yang bertahan hidup pascaperawatan. Namun, karena kami berfokus pada pola berulang dan bersama di antara pasien untuk mengidentifikasi jalur kemoresistom, kami yakin bahwa potensi bias apa pun karena heterogenitas tumor berdampak minimal pada hasil kami. Meskipun demikian, eksperimen sel tunggal di masa mendatang atau uji transkriptomik spasial akan sangat berharga untuk melengkapi dan memvalidasi temuan kami.

5 Kesimpulan
Singkatnya, penjelasan prinsip-prinsip yang mengatur mekanisme resistensi adaptif terhadap kemoterapi telah difasilitasi oleh kumpulan data longitudinal unik kami yang digabungkan dengan analisis klasifikasi pola gen. Adaptasi fenotipik yang kompleks dan heterogen yang digambarkan dari penelitian ini tampak sebagai kendala utama untuk pengobatan antikanker yang berhasil. Meskipun demikian, konsep-konsep mendalam yang diperoleh dari penelitian ini dapat mempercepat pengembangan metode yang tersedia secara klinis untuk memperkirakan keadaan adaptif yang muncul, yang pada akhirnya mengarah pada optimalisasi rencana pengobatan khusus pasien di masa mendatang.