השרדות: איך אפשר לדעת מה יהיה?

בעקבות הערתו של עופר לרשימה בה נתתי לינק למצגת ההרצאה שלי בכנס האחרון של האיגוד הישראלי לסטטיסטיקה, אני מפרסם כאן הסבר כללי ולא טכני שמיועד לקהל הרחב על הנושא בו עסקה ההרצאה. אני חייב להתריע בפני הקוראים כי מדובר ברשימה ארוכה מהרגיל, ולמרות כל מאמצי ייתכן מאוד שהיא עדיין לא מספיק ברורה עבור הקהל הרחב (אני מקווה מאוד שכן). בכל מקרה, אני שמח לפתוח כאן חלון אל העולם שבו אני פועל.

אחד הנושאים שמעניינים אותי הוא ניתוח השרדות, ותכנון של ניסויים קליניים לבדיקת אפקט של טיפול תרופתי על ההשרדות.

מה זה אומר? במחלות רבות, משך הזמן שעובר עד להתרחשותו של אירוע קליני מסויים (ומוגדר מאש) מהווה אינדיקציה לחומרת המחלה. במחלות כמו סרטן ו-ALS (הידועה גם בשם “מחלת לו גריג“, החולה המפורסם ביותר כיום במחלה זו הוא הפיזיקאי סטיבן הוקינג) האירוע הוא בדרך כלל מוות. בטרשת נפוצה, מחלה שמעניינת אותי במיוחד, אירוע כזה יכול להיות התקף נוירולוגי או הדרדרות ביכולות הפיזיות (על פי סולם ה-EDSS). מסיבות מובנות יש עניין לבדוק האם טיפול תרופתי (קיים או חדש) , יכול להאריך את משך הזמן הממוצע עד לאירוע, או, במילים אחרות, להקטין את קצב התרחשות האירועים.

איך בודקים את זה? אפשר כמובן לקחת קבוצת חולים, לתת לחלק מהם את התרופה שאנו רוצים לבחון ולחלק מהם טיפול ביקורת (פלסבו או הטיפול הסטנדרטי המקובל), ולהמתין עד שכל אחד מהחולים יעבור את האירוע שבו מתעניינים. אז נוכל להשוות את משך הזמן הממוצע עד להתרחשות האירוע בקבוצה אחת עם מקבילו מהקבוצה השניה. ניסוי כזה הוא בעייתי מאוד. הוא עלול להמשך זמן רב מאוד, ולמעשה משכו אינו מוגדר מראש. לכך יש השלכות מעשיות ואתיות כאחד. איזה חברת תרופות תתחיל ניסוי שלא ידוע מראש מה יהיה משכו ומה תהיה עלותו? איזה חולה יסכים להסתכן בטיפול בפלסבו לזמן בלתי מוגבל.

אפשר לתכנן את הניסוי באופן אחר: נמתין פרק זמן מסויים, קבוע מראש, נניח שנה. חלק מהחולים יעברו בפרק זמן זה את האירוע הקליני, ועבורם נדע בדיוק את משך הזמן שעבר מתחילת הטיפול ועד האירוע. שאר החולים לא יחוו את האירוע, ועבורם לא נדע את משך הזמן המדויק עד התרחשות האירוע, אבל כן נדע כי הוא ארוך משנה. מתברר שהאינפורמציה הזו מספיקה כדי לאמוד את משך הזמן הממוצע על לאירוע וכן לקבוע האם יש הבדלים בין שתי קבוצות הטיפול. יותר מכך, אפשר לקבוע מראש את גודל המדגם הדרוש כדי לזהות הבדל באופן מובהק סטטיסטית.

גודל המדגם הדרוש קשור ישירות למספר האירועים הצפויים להתרחש במהלך הניסוי, והקשר הוא דרך ההנחות לפיהן הניסוי מתוכנן. באחד השקפים הבאתי דוגמא שבה משך הזמן הממוצע עד לאירוע בחולים הלא מטופלים הוא בערך שש שנים, ורוצים לבדוק האם טיפול כלשהו מאריך את משך הזמן הממוצע ל-8 שנים בערך. חישבתי ומצאתי כי מדגם של 1000 חולים, שיטופלו במשך שנתיים, יספק עצמה סטטיסטית של 83% לזיהוי הבדל מובהק בין הקבוצות. אם מניחים בנוסף כי הסיכון לאירוע הוא קבוע לאורך הזמן, אפשר לחשב ולמצוא כי בתנאים האלה בערך 30% מהחולים הלא מטופלים (כלומר, מהחולים המטופלים בפלסבו) יעברו אירוע, כלומר 150 חולים בערך. לעומת זאת, רק 22% מהחולים המטופלים בטיפול הנסיוני אמורים לעבור אירוע (אם התרופה באמת עובדת), כלומר רק 110 מבין 500 החולים המטופלים יעברו את האירוע. לכן, בסופו של דבר אנו אמורים לצפות ב-260 אירועים. אפשר לעשות את החשבון הזה גם מהסוף להתחלה: תחילה לחשב כי יש צורך לצפות ב-260 אירועים כדי להגיע לעצמה הסטטיסטית הדרושה, ומשם לקבוע את גודל המדגם (1000 חולים).

המשמעות של כל הדיון האחרון: אם במהלך הניסוי ייצפה מספר קטן יותר של אירועים מכל סיבה שהיא, אזי העצמה הסטטיסטית של הניסוי תהיה נמוכה יותר. ולכן, השאלה שמדירה את שנת הסטטיסטיקאי (וכל צוות הפרויקט, בעצם) היא: האם בסוף הניסוי נגיע למספר האירועים הרצוי?

בואו נעבור לרגע ללוגיסטיקה של הניסוי הקליני. ברור לחלוטין כי אי אפשר לגייס בבת אחת 1000 חולים שיסכימו להשתתף בניסוי ויתחילו כולם את הטיפול באותו יום. גיוס של כמות כה גדולה של חולים יכול לארוך זמן רב. נניח שמשך הגיוס הוא שנה וחצי, כלומר בממוצע מצטרפים לניסוי בכל יום  כשני חולים ( 1.828 חולים, למתעניינים). פירוש הדבר הוא שאחרי שנה ורבע יהיו בתוך הניסוי כ-830 חולים. הראשון בהם כבר מטופל שנה ורבע, האחרון שהצטרף לניסוי התחיל את הטיפול זה עתה. מכאן מקבלים כי משך זמן הטיפול הממוצע לחולה בנקודת זמן זו הוא 0.625 שנים (כשבעה וחצי חודשים). חלק מהחולים האלה כבר עברו מן הסתם את האירוע הקליני בו אנו מתעניינים. נניח שעד עתה אירעו בסך הכל 60 אירועים. ובכן – האם בקצב הזה נגיע אל 260 האירועים הדרושים להצלחת הניסוי? איך נוכל לדעת זאת, כאשר איננו יודעים איזה טיפול מקבל כל אחד מהחולים?

בגלל השאלה הזו, גבירותיי ורבותיי, אני לא ישן טוב בלילה. אם אינכם יודעים, ניסויים קליניים משרים קונצנזוס. אם הניסוי מסתיים בהצלחה, אזי יש תמימות דעים בשאלה מי אחראי להצלחה. כל אחד חושב שהכל רק בזכותו. לעומת זאת, אם הניסוי מסתיים בכשלון, גם אז יש תמימות דעים: כולם מאשימים את הסטטיסטיקאי. ולכן, אחד מתפקידיו החשובים של הסטטיסטיקאי הוא לעקוב אחרי הנתונים המצטברים במהלך הניסוי ולהתריע על כל בעיה שהוא יכול לזהות.

אך איך אפשר לנסות להעריך כמה חולים ישרדו עד סוף הניסוי וכמה יעברו את האירוע?

יש תשובה נאיבית לשאלה הזו. כאמור, יש בניסוי כרגע 830 חולים עם משך טיפול ממוצע של 0.625 שנים, כלומר עד עתה צפינו במשך הניסוי ב-520 (830×0.625) שנות טיפול בערך. אנו מצפים בסך הכל לצפות ב-2000 שנות טיפול (1000×2), כלומר עד עכשיו צפינו ב-26% בערך מסך החשיפה הצפויה לטיפול. אם מספר האירועים פרופורציוני לסך החשיפה, הרי ש-60 האירועים שבהם צפינו עד כה מהווים 26% ממספר האירועים הצפוי, ולכן בסך הכל צפויים כ-230 אירועים בסוף הניסוי. לא טוב, אבל יש מה לעשות: אפשר להגדיל את המדגם ולגייס עוד 127 חולים מעבר ל-1000 שתכננו (תרגיל: איך הגעתי ל-127?). אבל האם ההנחה לפיה מספר האירועים פרופורציוני לחשיפה נכונה?

כאשר כתבתי את החלק הסטטיסטי של פרוטוקול אחד הניסויים הקליניים שהחברה בה אני עובד עורכת, התלבטתי רבות בשאלה הזו. הפתרון הנאיבי שתיארתי זה עתה, נראה לי, איך לומר, נאיבי. הקדשתי מחשבה לנושא ומצאתי פתרון שונה. התייעצות עם פרופ’ פול פייגין מהטכניון הביאה לפתרון שמבוסס על גישה אחרת, ובסופו של דבר זה הפתרון שנכלל בפרוטוקול.

חלף זמן מה, ושינוי ארגוני בחברה שינה את מבנה ההנהלה, וכך נתבקש סטטיסטיקאי העובד בשלוחת החברה בארה”ב, ד”ר יויי שו לחוות את דעתו על הפרוטוקול. אותו סטטיסטיקאי, הפנה אותי לשני מאמרים שפרסם הסטטיסטיקאי האנגלי ג’ון וויטהד, בהם גם הוא ניסה לענות על השאלה הזו. במאמר שפרסם וייטהד בשנת 2001 הוא הציע מתווה לחיזוי מספר האירועים הצפויים בסוף ניסוי כזה, על סמך הנתונים החלקיים שקיימים במהלך הניסוי: משך הזמן שבו כל חולה מטופל, האם החולה עבר אירוע, ואם כן, מתי, אך מבלי לדעת באיזה טיפול מטופל כל חולה (בהתאם לעקרון הסמיות הכפולה שעומד בבסיס כל ניסוי קליני). במתווה של וייטהד יש שלושה שלבים. תחילה יש לאמוד את קצב האירועים הכללי באוכלוסיית הניסוי, שהוא עירוב של קצב האירועים באוכלוסיית החולים המטופלים בתרופה הנסיונית וקצב האירועים באוכלוסיית חולי הפלסבו. בשלב השני יש להפריד באופן כלשהו את הקצב הכללי לשני הקצבים של שני הקבוצות, תוך הסתמכות על ההנחות ששימשו לתכנון הניסוי. לאחר ההפרדה, אפשר לחשב את מספר האירועים הצפוי מתוצאות השלב השני.

נשמע כמו תכנית טובה, והשלב הראשון והשלישי פשוטים למדי, אבל איך מיישמים את השלב השני? את זה וייטהד הסביר במאמר נוסף שפרסם בשנת 2001, עם שותפים. הפתרון שהוא הציע שונה משני הפתרונות ששקלתי בעת כתיבת הפרוטוקול. אבל החדשות המשמחות מבחינתי הן ששני הפתרונות ששקלתי משתלבים היטב במתווה של וייטהד, והם פתרונות קבילים על פי המאמר.

מאחר ויש בידנו שלושה יישומים אפשריים, עלתה השאלה האם יש עוד דרכים ליישם את המתווה של וייטהד? מסתבר שכן. התייעצות עם פרופ’ דוד צוקר מהאוניברסיטה העברית הניבה פתרון שהוא הרחבה של הפתרון המקורי שלי, והתייעצות עם ד”ר עפרה ברנט שעובדת יחד איתי הביאה לפתרון אפשרי נוסף. הדרך שהציע ד”ר שו התגלתה כשקולה לזו שמופיעה בפרוטוקול.

עכשיו, כשבידי 5 גישות אפשריות לאמידת מספר האירועים הצפוי בסוף הניסוי, השאלה הבאה שעלתה ברורה מאליה: איזו גישה היא הטובה ביותר? אני יכול, כמובן, לערוך את החישובים על פי 5 הגישות, ולהמתין לסוף הניסוי. דרך מהירה יותר היא לערוך סימולציה של ניסוי קליני, בתנאים שונים, ולראות איזה גישה חוזה את מספר האירועים באופן הטוב ביותר.

בחרתי לערוך סימולציה של הניסוי שתיארתי קודם: 1000 חולים, בשתי קבוצות טיפול, כאשר משך הטיפול הוא שנתיים, ומשך גיוס החולים הוא שנה וחצי, כך שמשך הניסוי ממועד גיוס החולה הראשון ועד סיום הטיפול של החולה האחרון הוא 3.5 שנים. הסימולציה קבעה לכל חולה את משך הזמן עד לאירוע באמצעות התפלגות וויבול (הסבר מייד). אם משך הזמן שהוגרל היה גדול משנתיים, החולה נחשב כשורד. כמו כן, חישבתי את הזמן הקלנדרי של האירוע, ואם הוא התרחש במשך השנה ורבע הראשונות של הניסוי, הוא נחשב כאירוע שצפיתי בו לצורך חיזוי מספר האירועים בסוף הניסוי.

מהי התפלגות וויבול? זוהי התפלגות סטטיסטית גמישה למדי, שמקובלת כמודל טוב עבור נתוני השרדות. להתפלגות זו שני פרמטרים. האחד מבטא את הסיכון של החולה לעבור אירוע קליני (מסומן באות היוונית ביתא במצגת). כאשר ביתא נתון, הפרמטר השני מאפשר לבטא את הסיכוי שחולה ישרוד עד נקודת זמן מסויימת (כלומר, לא יעבור אירוע קליני במשך תקופה מסויימת).

ערכו של הפרמטר ביתא משקף את אופי המחלה. תארו לכם רופא, חולה וסטטיסטיקאי שנפגשים יחד (ברכבת, או בבית המשפט, לא ממש חשוב איפה). החולה מספר כי הוא כבר X שנים במצבו הנוכחי, ולכן הוא אופטימי, מאמין כי מצבו ימשיך להישאר יציב, ובמלים אחרות, הוא מאמין כי ככל שחולף הזמן, הסיכון להחמרה במחלתו הולך וקטן. הרופא ממהר לצנן את התלהבותו, ואומר לא כי אם עבר הרבה זמן ללא הדרדרות, זהו אות לכך שהאירוע הקליני הולך וקרב. כלומר, הרופא מאמין כי הסיכון להחמרה במחלה הולך וגדל עם הזמן. הסטטיסטיקאי עובד בחברת תרופות, ראה אלפי נתונים של אלפי חולים, ולכן הוא מחייך לעצמו. הנתונים שראה מעידים, לדעתו, כי הסיכון הוא קבוע לאורך זמן, ומשך הזמן בו החולה נמצא במצב יציב לא מעלה ולא מוריד. (במחלת הטרשת הנפוצה זה כנראה נכון). ערכו של הפרמטר ביתא קובע את כיוון הסיכון: עולה (ביתא גדול מ-1), יורד (ביתא קטן מ-1) או קבוע לאורך זמן (ביתא שווה אחד). לכן ערכתי סימולציות למצבים בהם ערכי ביתא הם 0.5, 1 ו-2. הניסוי תוכנן בהנחה של סיכון קבוע לאורך זמן (ביתא=1).

הפרמטר השני, גאמא, מבטא את הקצב בו מתרחשים האירועים, ובעזרתו ניתן לבטא את אחוז החולים שישרדו (לא יעברו אירוע קליני)  במשך שתי שנות הטיפול. הניסוי תוכנן בהנחה ש-70% מהחולים המטופלים בפלסבו ישרדו, וערכתי סימולציות לשלושה מצבים: 70% שורדים, 60% שורדים (אוכלוסיה “פעילה” – פחות שורדים, יותר אירועים קליניים) וגם 80% שורדים (אוכלוסיה לא פעילה).

בסך הכל ערכתי לכן סימולציות של תשעה תרחישים אפשריים, שאחד מהם הוא התרחיש בו הנחות התכנון מתקיימות. הנה התוצאות של התרחיש הזה:

בציר האופקי מופיעות 6 נקודות: Actual הוא מספר האירועים בפועל בסוף הניסוי, Naïve הוא מספר האירועים החזוי על פי הגישה הנאיבית, וכן הלאה. (ML היא הגישה שציינתי בפרוטוקול הניסוי). בציר האנכי מופיע מספר האירועים החזוי בסוף הניסוי (תוצאות מ-500 הרצות של הסימולציה). הקו האפקי מסמן את 260 האירועים הדרושים. הגרפים האדומים הם דיאגרמות קופסא (box plots) שפותחו על ידי ג’ון טוקי (מספר 9 ברשימת 15 הסטטיסטיקאים הגדולים). הקצה העליון של כל מלבן מסמן את הרבעון העליון של הנתונים, הקצה התחתון מסמן את הרבעון התחתון, והקו שבאמצע המלבן מסמן את החציון. סימן ה-+ מסמן את הממוצע.

מהגרף הזה לומדים שני דברים. ראשית, כל השיטות נותנות תוצאות דומות,בממוצע. שנית, אם ההנחות לפיהן הניסוי תוכנן אכן נכונות, אז אין מה לדאוג.

אולם שמונת התרחישים האחרים פחות מרגיעים (מי שמעוניין להתעמק בהם מופנה אל המצגת). עדיין נראה כי כל חמשת השיטות נותנות תוצאות דומות בדרך כלל (עם כמה חריגות, שעדיין לא ברורה לי משמעותן). זה אומר שהשיטה הנאיבית היא לא כל כך נאיבית, או שמא השיטות המתוחכמות גם הן נאיביות למדי. אמנם, אם ההנחה של סיכון קבוע תקפה, אז השיטות חוזות באופן לא רע את מספר האירועים בניסוי (הניתוח שערכתי עד כה הוא פשטני ויש מקום לניתוח יותר מעמיק שבכוונתי לערוך), אבל במקרים האחרים הניבויים סוטים מהתוצאה האמיתית וחושפים אותנו לסיכונים. כרגיל בסטטיסטיקה יש שתי טעויות אפשריות. טעות אפשרית אחת היא ניבוי חסר של מספר האירועים הצפוי. עקב כך, ייתכן כי יוחלט להגדיל את המדגם (כלומר לגייס עוד חולים לניסוי, ולהאריך את משך הניסוי ועלותו) ללא צורך. ניבוי עודף של מספר האירועים עלול להביא להמשך הניסוי ללא שינוי, ורק בסופו יתברר כי לא התרחשו מספיק אירועים והניסוי נכשל עקב כך.

מה הלאה? לאחר ההרצאה שאל אותי פרופ’ דוד שטיינברג מאוניברסיטת תל-אביב האם ניתן לבדוק אם הנחות הניסוי תקפות במהלך הניסוי. עניתי כי אני חושב שלא, אבל אחר כך עלה בדעתי רעיון שאולי יאפשר בדיקה האם הסיכון עולה, יורד או קבוע במהלך הניסוי, וזאת על ידי הרחבת הגישות שבדקתי לאמידה של מספר פרמטרים (במלים אחרות, אני רוצה לאמוד בו זמנית גם את מספר האירועים הצפוי וגם את הפרמטר ביתא). תיאורטית זה לא מסובך. טכנית זה ידרוש מקסימיזציה של פונקציה לא לינארית של שני משתנים, דבר שדורש קצת עבודה. מקווה שאוכל לעשות אותה מתישהו.