מידול מתמטי של התפשטות מגפות

מבט נוסף על מקדם ההדבקה ואיך מחשבים אותו

פרופ׳ ניר גביש, הפקולטה למתמטיקה, טכניון

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

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

בסמסטר אביב יינתן קורס ״מודלים מתמטיים בחקר מגיפות״ (198008).

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

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

Graph of reproduction number — גרף של מקדם ההדבקה, מתוך דו״ח התפשטות יומי 31.1 של מרכז המידע ומרכז המידע והידע הלאומי למערכה בקורונה

פוסט זה יעסוק במקדם ההדבקה: איך הוא מוגדר, מה הוא אומר ואיך מחשבים אותו?

מקדם ההדבקה הבסיסי ומקדם ההדבקה האפקטיבי

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

עם התפתחות המגפה, הסביבה של חולה ממוצע מורכבת מיותר ויותר נדבקים או מחלימים, ולכן האפשרות שלו להדביק אחרים קטנה. בהתאם, מקדם ההדבקה קטן עם התפתחות המגפה ועם הצטברות מחלימים. מקדם ההדבקה העדכני שלוקח בחשבון את הצטברות המחלימים והמחוסנים נקרא מקדם ההדבקה האפקטיבי. נהוג לסמן את מקדם ההדבקה האפקטיבי ב $R, R_t$ או $R_{\rm eff}$. במודל של קבוצת אוכלוסייה אחת מקדם ההדבקה האקטיבי נתון ע״י הנוסחא

$R_{\rm eff}=R_0S(t)$

כאשר $S(t)$ הוא מספר המועדים היחסי בזמן t (ראו פוסט לפרטים). בזמן אפס, כמעט כל האוכלוסייה מועדת ולכן $Sֿ(0)\simeq 1$ ו $R_{\rm eff}\simeq R_0$.

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

חישוב מקדם ההדבקה האפקטיבי

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

המשוואה הבאה מתארת את השינוי במספר הנדבקים הכולל

$\underbrace{I^\prime(t)}_{\large\rm שינוי\, במספר\, הנדבקים\, הכולל}=\underbrace{\frac1d R_{\rm eff}(t) I(t)}_{\large\rm מספר \,נדבקים\, חדשים\, ביום}-\underbrace{\frac1dI(t)}_{\large\rm מספר\, מחלימים\, חדשים \,ביום}$ (*)

כאשר $I(t)$ הוא מספר הנדבקים (החולים) הפעילים, $S(t)$ הוא מספר המועדים, $R_{\rm eff}(t)$ הוא מקדם ההדבקה האפקטיבי ו $d$ הוא מספר הימים עד להחלמה.

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

$R_{\rm eff}(t)=\frac{{\Large\{\rm מספר\, נדבקים\, חדשים\}}}{\frac1dI(t)}$

המכנה מתאר את סך החולים הפעילים חלקי משך המחלה בימים $d$.

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

$(*)\qquad R_{\rm eff}(t)=\frac{{\Large\{\rm מספר\, נדבקים\, חדשים\}}}{{\Large\{\rm פוטנציאל\, ההדבקה\, של\, החולים\, הפעילים\}}}$

לדוגמא, אם היום פוטנציאל ההדבקה להיום שקול ל 100 חולים פעילים ומהם נדבקו היום 120 אנשים, אז מקדם ההדבקה הוא 1.2. כל חולה פעיל הדביק בממוצע 1.2 איש.

החישוב בפועל

הנוסחא (*) משמשת לחישוב מקדם ההדבקה האפקטיבי, אך על מנת ליישמה בפועל יש להתחשב בשני גורמים עיקריים:

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

נתחיל מטיפול מדויק יותר בפוטנציאל ההדבקה של החולים הפעילים. אדם יכול להידבק מחולה שנמצא בשלבים שונים של המחלה שלו, נסמן ב $I_{t-1}$ את החולים שנדבקו לפני יום. באופן דומה, נסמן ב $I_{t-2}$ את החולים שנדבקו לפני יומיים וב $I_{t-s}$ את החולים שנדבקו לפני $s$ ימים. פוטנציאל ההדבקה הוא

אחוז ההדבקות שמתרחשות ע״י חולים יום אחרי ההידבקות שלהם מסומן ב $w_1$. באותו אופן $w_2$ הוא אחוז ההדבקות שמתרחשות ע״י חולים יומיים אחרי ההידבקות שלהם ו $w_s$ הוא אחוז ההדבקות שמתרחשות $s$ ימים לאחר ההידבקות. נהוג לחשוב על הנדבקים החדשים כדור נדבקים חדש, ילדים ונכדים לדורות הנדבקים הקודמים. בהתאם, הפונקציה $w_s$ המתקבלת נקראת התפלגות זמן הדור ( generation time distribution או gtd) ומגדירה את התפלגות הזמן בין המועד בו אדם נדבק לבין המועד בו הדביק אדם אחר. הצבת הביטוי לפונטציאל ההדבקה בנוסחא (*) נותן את הנוסחא הבאה למקדם ההדבקה (סימן ה $ֿ\Sigma$ מתאר סכימה עם ערכי $s$ עולים).

הנוסחא לחישוב R — מפרסומי מרכז המידע והידע הלאומי למערכה בקורונה

נעבור כעת לטפל באי-ודאות לגבי מספר הנדבקים. היות ומספר הנדבקים מופיע גם במונה וגם במכנה, נוכל להציב בנוסחא כל ביטוי שהוא יחסי לו. לדוגמא, אם בכל יום מאתרים מחצית ממספר הנדבקים, אז נוכל להציב בנוסחא עבור $R$ את מספר המאומתים החדשים ביום $D_t$

$R=\frac{I_t}{\sum_s I_{t-s}\times w_s}=\frac{D_t}{\sum_s D_{t-s}\times w_s}$

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

$R=\frac{I_t}{\sum_s I_{t-s}\times w_s}=\frac{D_{t+\tau}}{\sum_s D_{t+\tau-s}\times w_s}$

כאשר $\tau$ הוא זמן ההשהייה, כשבעה ימים. באופן דומה ועם זמן השהייה מותאם, ניתן לחשב את מקדם ההדבקה בעזרת נתוני חולים או נתוני תמותה.

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

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

מקדם ההדבקה הנמדד מול מקדם ההדבקה בפועל

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

חיזוי לעכשיו וחיזוי לאחור

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

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

הפניות ומידע נוסף

מסמכי מרכז מידע והידע הלאומי למערכה בקורונה:

הקוד למודל החיזוי פתוח וזמין כאן.

בסמסטר אביב יינתן קורס ״מודלים מתמטיים בחקר מגיפות״ (198008).

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

זה המקום עבורכם להגיב ולשאול

סיום המגפה ולמה חיסוניות עדר לא תגן על מי שלא יתחסן

פרופ׳ ניר גביש, הפקולטה למתמטיקה, טכניון

אני מודה לפרופ׳ חגי כתריאל על הערות מועילות.

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

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

מהלך המגפה וסף חיסוניות העדר

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

נתייחס לשני תרחישי סיום שונים. שניהם מתחילים בסגר הדוק, היות והסגר השלישי הוא עובדה מוגמרת בשלב זה.

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

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

תרחיש א׳ מתאר חציית סף חיסוניות העדר ע״י חיסון מסיבי של האוכלוסייה. זו דוגמא קלאסית להגנה באמצעות חיסוניות עדר על אלו שאינם יכולים או אינם רוצים להתחסן. האם התרחיש ריאלי? הנוסחה לסף חיסוניות העדר* היא $1-1/R$ (מפנה שוב לפוסט הזה לפרטים), כאשר עבור R=2.5 המיוחס לזן הקורונה המוכר סף חיסוניות העדר הוא 60%. עם זאת, הווריאנט הבריטי (והדרום-אפריקאי) נמצא כעת בתהליך חדירה וסביר מאוד להניח שיהיה הזן הדומיננטי בשלבי הסיום של המגפה. מקדם ההדבקה של הווריאנט הוא סביב ה 4 ובהתאם סף חיסוניות העדר עבורו הוא סביב 75%. כ 30% מאוכלוסיית ישראל נמצאת מתחת לגיל 16 ולכן לא יכולה להתחסן, בנוסף ניתן להעריך שכ 10% מהאוכלוסייה כבר חלתה ופיתחה נוגדנים לנגיף. המשמעות היא שיש לחסן את כל האוכלוסייה הבוגרת כדי להגיע לסף חיסוניות עדר או קרוב אליה. זה יעד היענות שהוא בלתי מושג. במידה ו 80% מהאוכלוסייה הבוגרת תתחסן, סף חיסוניות העדר ייחצה רק כאשר כ 15% נוספים מהאוכלוסייה הכללית הלא-מחוסנת (בוגרת וצעירה) תידבק. ללא ניהול נכון של התהליך, מדובר ביותר מהכפלה של התחלואה המצטברת והתמותה עד כה, וזאת על סמך ההנחה שהתחלואה והתמותה עד כה נבעה מהידבקות של 10% מהאוכלוסייה. מספרים אלו ישתנו לרעה במידה והחיסון יעיל רק חלקית במניעת הדבקה.

תרחיש ב׳ מתאר למעשה התפרצות מבוקרת של המחלה בה התחלואה הקשה והתמותה נשלטות ע״י חיסון האוכלוסייה הבוגרת ובפרט אוכלוסיית הסיכון. לאחר חציית סף חיסוניות העדר, המגפה דועכת באיטיות ובזמן הדעיכה עוד ועוד אנשים נדבקים והתחלואה מצטברת. מה אחוז האוכלוסייה שתידבק בסופו של דבר? תהליך הדעיכה מתואר ע״י מודל SIR (מפנה שוב לפוסט לפרטים) והנוסחא לאחוז הנדבקים הכולל* $S_\infty$ היא $$ S_\infty=1+\frac{1}{R_{0}}\log S_\infty $$ עבור R=4 המיוחס לווריאנט הבריטי, 98% מהאוכלוסייה תידבק. יעילות חלקית (אפילו 95%) של החיסון במניעת הדבקות תעלה את המספר הזה כלפי מעלה. המשמעות היא שיש לקחת בחשבון שכל מי שלא מתחסן, יידבק.

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

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

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

סימולציה של התפתחות תחלואה עם 83% היענות לחיסון מקרב בני 60 ומעלה

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

סימולציה של התפתחות תחלואה עם 83% היענות לחיסון מקרב בני 60 ומעלה

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

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

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

מצגת לצט״מ (צוות טיפול מגפות) בנושא זמינה כאן

הקוד למודל התחלואה פתוח וזמין כאן.

זה המקום עבורכם להגיב ולשאול

הם הסגר יכול להאיץ את חדירת הווריאנט הבריטי?

האם הסגר יכול להאיץ את חדירת הווריאנט הבריטי?

פרופ׳ ניר גביש, הפקולטה למתמטיקה, טכניון

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

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

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

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

נתחיל מתהליך ההתפשטות וההתבססות של הווריאנט הבריטי ללא סגר. בכל החישובים הבאים נניח שאופי המחלה – תקופת דגירה, זמן להחלמה וכן הלאה – זהה בשני הזנים, אך בווריאנט הבריטי חולה מדביק בממוצע 50% יותר אנשים במשך תקופת המחלה שלו. המשמעות היא שלכל זן יש מקדם הדבקה משלו. עבור הזן המוכר נתייחס למקדם הדבקה של 1.2 ועבור הווריאנט הבריטי 1.8 (50% יותר). אלו הם מקדמי הדבקה שלוקחים בחשבון הגבלות שונות ושמירה על ריחוק חברתי, פחות או יותר בטווח מקדם ההדבקה הנוכחי (טרם הסגר ההדוק). אם בזמן אפס יש 1000 נדבקים בזן המוכר ו-10 נדבקים בווריאנט הבריטי אז כל חולה בזן המוכר ידביק 1.2 אנשים ונקבל 1200 נדבקים חדשים בזן המוכר עד לסיום תקופת המחלה שלהם. במקביל, כל חולה עם הווריאנט הבריטי ידביק 1.8 אנשים באותה תקופה ונקבל 18 נדבקים חדשים. בשבוע העוקב, 1200 החולים החדשים ידביקו 1444 אנשים בזן המוכר ואותם 18 חולים עם הזן הבריטי ידביקו 32 אנשים. כך משבוע לשבוע, מספרי הנדבקים יעלו ואחוז הנדבקים עם הווריאנט הבריטי יטפס עד שיהפוך לזן הדומיננטי.

הנוסחא תיראה כמו

$ (*) ֿ\qquad I_1(t)=I_1(0) R0^{t/\tau},\quad I_2(t)=I_2(0) (1.5R0)^{t/\tau}$

כאשר $ \tau$ הוא פרק הזמן מהידבקות של אדם ועד להדבקת אחרים (סדר גודל של 5 ימים בקורונה). הטבלה הבאה מסכמת את נתוני התפשטות לאורך מספר תקופות זמן

זמן	0	$ \tau$	$ 2\tau$	…	$ 6\tau$	$ 10\tau$	$ 15\tau$
נדבקים בזן המוכר $ I_1(t)$	1000	1200	1444		3000	6200	15,400
נדבקים בווריאנט הבריטי $ I_2(t)$	10	18	32		340	3600	67,500
אחוז נדבקים בווריאנט הבריטי $ \frac{I_2(t)}{I_1(t)+I_2(t)}$	0.9%	1.5%	2.2%		10%	36%	81%

מה קורה בזמן סגר? נניח שבסגר מקדם ההדבקה של הזן המוכר יורד ל 0.8. מקדם ההדבקה של הווריאנט הבריטי גבוה ב 50% ולכן הוא 1.2 בזמן אותו הסגר. המשמעות היא שמספר הנדבקים בזן המוכר קטן בזמן הסגר, בעוד שמספר הנדבקים בווריאנט הבריטי גדל. האם במצב כזה הסגר מאיץ את חדירת הווריאנט הבריטי? לא. מנוסחא (*) רואים שאחוז הנדבקים בווריאנט הבריטי אינו תלוי ב R

$ \frac{I_2(t)}{I_1(t)+I_2(t)}=\frac{1.5^{t/\tau} I_2(0)}{I_1(0)+1.5^{t/\tau} I_2(0)}$

לכן קצב החדירה אינו מושפע מהשינוי ב R בעקבות הסגר. הסגר לא יאט את חדירת הווריאנט הבריטי ולא יאיץ אותו.

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

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

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

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

זה המקום עבורכם להגיב ולשאול

האם החיסונים יכולים לעצור את גל התחלואה הנוכחי?

האם חיסונים יכולים לבלום את גל התחלואה הנוכחי?

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

חישוב זריז על גב מעטפה

עד כה הצגתי בבלוג משוואות דיפרנציאליות וחישובים אחרים שנעשים במחשב, אבל דווקא החישובים הזריזים, חישובי back of the envelope, נותנים תמונה ראשונית ומבססים הבנה. הנה גרסה של חישוב נפוץ להשפעת מבצע החיסונים. כ 70% מהחולים הקשים בארץ מגיעים מקבוצה של כ 1.4 מיליון בני השישים ומעלה החיים בארץ. בקצב חיסון של 50,000 מחוסנים חדשים ביום (שווה ערך, יחד עם המנה השנייה, ל 100,000 חיסונים ביום החל מהשבוע השלישי למבצע החיסונים), ניתן לחסן כרבע מקבוצה זו בשבוע. לכן בהשהיה מסוימת, כל שבוע חיסונים יפחית 17.5% מהתחלואה הקשה בארץ, יפחית עומס ממערכת הבריאות ואולי ייתר את הסגר המתממש. תוך חודשים ספורים, לאחר שמבצע החיסונים יסתיים והחיסון יגיע ליעילות גבוהה, הקורונה תהיה מאחורינו, לפחות כאירוע שעלול לצאת משליטה. נחזור לחישוב הזה ולמסקנות בסוף הפוסט.

תרגיל חימום

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

לסיכום, ניתן לראות מהגרף:

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

תרגיל החימום הוא בעצמו חישוב או טיעון זריז שנותן הבנה ראשונית של השפעת החיסונים על התחלואה. עם זאת, הוא מוגבל בכך שאינו אומר הרבה על צורת הגרף בטווח $1<R<1.3$. בפרט, הוא נעדר פרטים על יעילות החיסון, שינוי ביעילות לאחר המנה השנייה, ואופן הקצאת החיסונים. הטיעון גם לא יכול להתייחס למצב של שינוי במגבלות או התנהגות הציבור. ניתן כעת להציע חישוב מעט יותר מורכב שלוקח חלק מהפרטים האלו בחשבון, אך במקום זה אני רוצה לקפוץ ישירות לחיזוי של מודל מלא.

חישוב מלא

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

ניתן לראות שבקצב של 50,000 מתחסנים חדשים ביום ניתן לבלום את גל התחלואה ע״י מבצע חיסונים כל עוד מקדם הדבקה לא גדול מ 1.15. כצפוי יש רגישות גבוהה למקדם ההדבקה – במקדם הדבקה של 1.2 נדרשים לקצב של 100,000 מתחסנים חדשים ביום (יחד עם המנה השנייה – 200,000 חיסונים ביום) ובמקדם הדבקה של 1.25 נדרשים לחיסון של 300,000 מתחסנים חדשים ביום.

מסקנות וחזרה לחישוב הזריז מתחילת הפוסט

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

מה לגבי השאלות הקונקרטיות לגבי החיסונים וגל התחלואה הגואה? במקדם הדבקה 1.25-R=1.3 שבו אנחנו נמצאים, למיטב הבנתנו וע״פ תחזיות המודל, גל התחלואה הנוכחי לא ייבלם רק בזכות מבצע החיסונים, אפילו אם נגדיל מאוד את קצב החיסון. סעיף 3 של סיכום קבינט המומחים מפרט תרחישים שונים ודרכי פעולה עיקריות שנבחנו במודל החיסונים, יחד עם המסקנות שעלו ממנו.

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

זה המקום עבורכם להגיב ולשאול

שכבות גיל ומטריצות אינטרקציות

שכבות גיל ומטריצות אינטרקציות

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

הקדמה – חלוקה לקבוצות

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

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

חלוקה לקבוצות גיל

פוסט זה יתמקד בהרחבה לקבוצות גיל והצגת מטריצות האינטראקציות שצצות בהרחבה כזו. מודלים רבים מתייחסים ל 16 קבוצות גיל – קבוצת גילאים 0-4, 5-9, 10-14, וכן הלאה עד קבוצת גילאים 70- 74 וקבוצה אחרונה 75+. באופן כללי אפשר לדבר על $n$ קבוצות. אם במודל SIR בסיסי, המשתנים $S,I,R$ תיארו את כלל המועדים, הנדבקים והמחלימים באוכלוסייה, בהתאמה, עכשיו לכל קבוצה בנפרד יהיו משתנים $S_{i},I_{i},R_{i}$ שמתארים את מספר המועדים, הנדבקים והמחלימים בקבוצה $i$.

נסתכל עכשיו על קצב ההדבקה. ניזכר שבמודל SIR הבסיסי השינוי במספר המועדים שווה למינוס קצב ההדבקה וניתן ע״י הביטוי

$$S'(t)=-\beta S(t)I(t) $$

כאשר $\beta$ הוא מספר האנשים הממוצע שנדבק מדביק ביום כאשר הוא מוקף רק במועדים. הפרמטר $\beta$ עצמו תלוי באנשים בממוצע נמצאים במגע עם הנדבק ומה הסיכוי שמגע יוביל להדבקה. ניתן לרשום אותו בצורה

$$\beta=\beta_s c$$

כאשר $c$ הוא מספר המגעים הממוצע ו $\beta_s$ הוא הסיכוי להדבקה במגע בודד.

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

$$S_i'(t)=-\beta_s(c_{i1} I_1+c_{i2}I_2+…c_{in} I_n)S_i $$

כאשר $c_{ij}$ הוא מספר המגעים הממוצע בין אדם בקבוצה $i$ לאדם בקבוצה $j$, ותחת ההנחה שהסיכוי להדבקה במגע בין אנשים אחיד לכל הגילאים.
באופן דומה, נקבל ביטוי לקצב ההדבקה לכל הקבוצות $i=1,2,\cdots,n$

$$\begin{eqnarray}S_1′(t)&=-\beta_s(c_{11} I_1+c_{12}I_2+…c_{1n} I_n)S_1\\
S_2′(t)&=-\beta_s(c_{21} I_1+c_{22}I_2+…c_{2n} I_n)S_2\\
&\vdots\\
S_n'(t)&=-\beta_s(c_{n1} I_1+c_{n2}I_2+…c_{nn} I_n)S_n
\end{eqnarray} $$

קיבלנו אם כן $n^2$ איברים מהצורה $c_{ij}$ שמתארים את מספר המגעים הממוצע בין אדם מקבוצה $i$ לאדם בקבוצה $j$, שניתן לסדר בטבלה או מטריצה. זוהי מטריצת האינטרקציות.

הערת ביניים: מחקרים חדשים תומכים בכך שהסיכוי של ילד להידבק נמוך מהסיכוי של מבוגר להידבק, כלומר הסיכוי להידבק תלוי גיל $\beta_s=\beta_s(i)$. באופן דומה, גם הסיכוי של ילד להדביק נמוך מהסיכוי של מבוגר להדביק, כלומר יש לתקן את הביטוי כך ש $\beta_s=\beta_s(i,j)$ ותתקבל מטריצה נוספת. הרחבה זו יכולה להיות משמעותית בתיאור תרומת מערכת החינוך להתפשטות המגפה. עם זאת, לטובת פשטות ההצגה, לא אעסוק בה. לפרטים נוספים לגבי מחקר הדבקת והידבקות ילדים, ראו שקפי הרצאה של פרופ׳ יאיר גולדברג, טכניון, שעסק במחקר הזה.

איך מאכלסים את מטריצות האינטרקציות?

עד כאן הצגתי הרחבה של מודל ה SIR בה מסתכלים על כל האוכלוסייה כמקשה אחת למודל בו מסתכלים על קבוצות הגיל השונות באוכלוסייה. בדרך נוספו עוד פרמטרים שמתארים את האינטרקציה בין קבוצות הגיל השונות. בתיאור מקובל של 16 קבוצות גיל, נוספו 256 פרמטרים $c_{ij}$ שמסודרים בצורת טבלה או מטריצה. כמו בהרבה מקרים, קל לרשום את המשוואות, קשה יותר לאכלס את טבלאות הנתונים. צריך עכשיו לעסוק בשאלות כמו מה מספר האינטרקציות הממוצע בין נער בן 15 לאדם בן 38. איך עושים את זה?

ב-2004 , כנראה בהשפעת התפרצות ה SARS, האיחוד האירופאי החל במימון פרוייקט מחקר בשם POLYMOD שמטרתו (בתרגום חופשי) ׳שיפור מדיניות בריאות הציבור באירופה באמצעות מידול והערכה כלכלית של התערבויות לבקרת מחלות זיהומיות׳. נשמע כמו משהו שעשוי להיות שימושי כעת, לא? הפרוייקט החל לבחון מודלים של SIR עם קבוצות גיל ונתקל באותו קושי של אכלוס מטריצת האינטרקציות. הפרוייקט שם לו ליעד לאסוף נתונים ולאכלס את הטבלאות האלו, ואכן בעזרת סדר גודל של 40,000 שאלונים שמולאו בשמונה ארצות אירופאיות נותחו המגעים הבין-אישיים ונבנו מטריצות אינטרקציות. המידע שנאסף במחקר הזה הוא הבסיס, במקרים רבים, למטריצות אינטרקציות שמשתמשים בהם במודלים ברחבי העולם.

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

Projection of contact matrix — תהליך הטלה של נתונים מ POLYMOD למטריצת אינטרקציה ספיציפית למדינה. תמונה הותאמה מ https://doi.org/ 10.1371/journal.pcbi.1005697.

איך נראית מטריצת אינטרקציות?

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

ניתן לראות בגרפים המצורפים דוגמאות טיפוסיות של מטריצות אינטרקציה, כאשר צבע כהה מסמל מספר מגעים גדול וצבע לבן אפס מגעים. המטריצה (a) מתארת מגעים שמקורם בבית. ניתן לראות אלכסון ראשי דומיננטי שמשמעו נטייה של אנשים לגור או לארח אנשים בגיל דומה להם. שני האלכסונים המשניים בפער של כ 25-30 שנה נובעים מאינטרקציה בבית עם הילדים או עם הסבים והסבתות. המטריצה (d) מתארת מגעים שמקורם בעבודה. ניתן לראות שכעת האינטרקציות מוגבלות רק לגילאי העבודה 20-65. האינטרקציות די שטוחות בטווח הגילאים הזה, עם דומיננטיות חלשה של האלכסון בשל נטייה של אנשים בגיל דומה להתחבר, לדוגמא, בזמן הפסקות. מטריצת המגעים (g) שמקורם בבי״ס שונה מאוד – כאן יש אלכסון בגילאי בי״ס עם מעט אינטרקצייה מחוץ לשכבת הגיל שנובעות בעיקר ממגע עם מורים. לבסוף, המטריצה (j) מתארת מגעים שמקורם אינו בבית, בעבודה ובבי״ס. לדוגמא, בפעילויות ספורט, בתי כנסת, מסעדות וכו׳. המטריצה הזו שטוחה יחסית שמשמעה פיזור של המגעים בין כל שכבות הגיל.

Interaction matrices — דוגמאות למטריצות אינטרקציה למגעים ממקורות שונים. תמונה הותאמה מ https://doi.org/ 10.1371/journal.pcbi.1005697.

בחינת צעדי מדיניות תוך שימוש במטריצות האינטרקציות

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

לבסוף, באמצעות מטריצות האינטרקציות ניתן להעריך את תרומת צעדי ההתערבות לריסון המגפה. ניתן, לדוגמא, לחשב את מקדם ההדבקה במערכת עם מטריצת אינטרקציות שהופחתו ממנה מגעים שמקורם בבתי ספר ובכך להעריך את תרומת מערכת החינוך להתפשטות המגפה. בבדיקה ראשונית כזו, חישוב שלי כזה הראה שסגירת מערכת החינוך מפחיתה בכ 8% את מקדם ההדבקה. עם שקלול מחקר גרטנר לגבי הדבקת והידבקות ילדים, התרומה צנחה ל 4%-3% בדומה לנתונים אמפיריים שנאספו ב 11 ארצות באירופה.

אילו עוד רכיבים נכנסים למודל ומה עוד חסר?

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

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

רגע, זה לא המודל ששגה בגדול?

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

רגע, זה לא המודל ששגה בגדול?

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

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

Graph of simulator output vs data

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

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

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

האם התחזיות הראשוניות של עשרות אלפי מתים היו סבירות? כן ולא. לקחו מודל שמתאר בצורה טובה התפשטות של שפעת עונתית והתאימו אותו למאפיינים של קורונה. בשפעת עונתית, 20% מהאוכלוסייה נדבקת בשפעת, כלומר יש כ 1.8 מיליון נדבקים בישראל. בשיא העונה בישראל כ 8,000 איש ביום פונים לקופ״ח או למיון עם סימפטומים של שפעת והתפוסות במח׳ פנימיות בבתי״ח מגיעות לשיא. כל זאת כאשר כ 60% מהאוכלוסייה מעל גיל 60 מתחסנת (וכ 20% מהאוכלוסייה הכללית). אלו נתונים מבוססים שחוזרים על עצמם שנה אחר שנה (לפעמים בגל תחלואה קשה יותר ולפעמים פחות). ניתוח מאות מקרים מההתפרצות בסין הראה שהקורונה מדבקת פי כמה משפעת עונתית וקטלנית פי כמה ממנה. תמונה זו לא השתנתה מהותית לאחר מיליוני מקרים בעולם. מה קורה כאשר מזינים לאותו מודל של שפעת עונתית פרמטרים של שפעת ׳על סטרואידים׳? מקבלים תחזית של מיליוני נדבקים, עשרות אלפים חולים ביום שמבקשים טיפול, צרכי אשפוז מרקיעי שחקים וכן הלאה. צריך להתייחס מאוד בחשדנות למודל שנותן תחזית שונה.

Snapshot from Twitter regarding covid19 — השוואה של קורונה מול שפעת. סוף פברואר. טוויטר של אורי שליט, תעו״נ.

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

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

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

Models developed for COVID

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

זה המקום עבורכם להגיב ולשאול

מה ההשפעה של אי-ציות? מה ההשפעה של בידוד סלקטיבי?

מה ההשפעה של אי-ציות ושל בידוד סלקטיבי?

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

הקדמה (ארוכה מהרגיל)

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

מה קורה כאשר הציות להנחיות הוא חלקי? האם אי-ציות, נניח של $20\%$ מהאוכלוסייה, יוריד לטמיון את תוצאות מאמצי שאר האוכלוסייה לשמור על ההנחיות או שיתבטא רק בפגיעה קטנה?

אפשר גם לתהות מה קורה בכיוון ההפוך –

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

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

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

מודל SIR עם קבוצות

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

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

נסמן ב $S_{1},I_{1},R_{1}$ את הגודל היחסי של מחלקות המועדים/נשאים/מחוסנים של קבוצה 1 וב $S_{2},I_{2},R_{2}$ את הגודל היחסי של מחלקות המועדים/נשאים/מחוסנים של קבוצה 2. אסמן ב $\beta_{ij}$ את מספר האנשים הממוצע מקבוצה $i$ שנשא מקבוצה $j$ מדביק ביום. קצב ההדבקה שווה למספר המועדים הממוצע שנשא מקבוצה 1 מדביק ביום ועוד מספר המועדים הממוצע שנשא מקבוצה 2 מדביק ביום. המשוואה המתקבלת היא \[ S_{i}^{\prime}(t)=-\left[\beta_{i1}I_{1}(t)+\beta_{i2}I_{2}(t)\right]S_{i}(t),\qquad i=1,2, \]

(זה מאוד דומה לפיתוח המודל הבסיסי, לאלו שנעים בחוסר נוחות עכשיו אני ממליץ לחזור לפוסט הראשון שמתאר את פיתוח המודל הבסיסי)

בדומה למודל הבסיסי, קצב ההחלמה הוא \[ R_{i}^{\prime}(t)=\frac{1}{d_{i}}I_{i}(t),\qquad i=1,2, \]

ע”י סכימת המשוואות ניתן לקבל את המשוואה לשינוי במספר הנשאים בכל אחת מהקבוצות \[ I_{i}^{\prime}(t)=\left[\beta_{i1}I_{1}(t)+\beta_{i2}I_{2}(t)\right]S_{i}(t)-\frac{1}{d_{i}}I_{i}(t),\qquad i=1,2, \]

קיבלנו סה”כ 6 משוואות לשישה נעלמים, כאשר הגודל היחסי של הקבוצות באוכלוסייה אינו משתנה ולכן נקבע ע”י המצב ההתחלתי.

הצצה ראשונה לתוצאות

בואו נפתור את המשוואות באופן נומרי במחשב בתרחישים שונים ונקבל תשובות ראשונות לשאלות.

אתחיל משאלת הציות: קבוצה 1 תהיה קבוצת הצייתנים, היא נשמעת להנחיות ולכן נדבקת פחות ומדביקה פחות. הקבוצה השנייה תהיה קבוצת המזלזלים שלא נשמעת באותה קפדנות להנחיות ובהתאם בעלת סיכוי גבוה יותר להידבק ולהדביק. בהתאם נבחר, $\beta_{11}=0.1$ ו $\beta_{22}=0.3$, כלומר הסיכוי להידבקות במפגש בין נשא ומועד צייתנים ($\beta_{11}$) קטן פי שלושה מסיכויי הידבקות במפגש בין נשא ומועד מזלזלים ($\beta_{22}$). מה קורה במפגש בין נשא מזלזל ומועד צייתן? נבחר ערך ביניים $\beta_{12}=0.2$ וכנ”ל לגבי מפגש בין נשא צייתן ומועד מזלזל – כלומר $\beta_{21}=0.2$. חשוב לי להדגיש שאין נסיון להתאים את הערכים לנתונים, המטרה היא לענות על שאלת הציות, לאו דווקא עבור ערכי פרמטרים ספציפיים. אני מתחיל מבחירת הערכים כי זה נדרש להרצה במחשב, אבל אם הייתי מתחיל מאנליזה הייתי מדלג על רמת הפירוט בשלב הזה. הגרף הבא מציג את הגודל היחסי מכלל האוכלוסייה של מחלקת הנדבקים $I=I_{1}+I_{2}$ במקרה של אי-ציות חלקי של $5\%,20\%$ מול תרחיש של ציות מלא, כלומר כאשר קבוצה 2 היא בגודל יחסי של $0,0.05,0.2$, בהתאמה.

Simulation of model with two groups ניתן לראות שמספר הנדבקים בשיא קופץ ב $66\%$ במקרה של אי-ציות חלקי של $5\%$ מהאוכלוסייה ופי ארבעה במקרה של $20\%$ אי-ציות! מבחינה פרקטית, זה יכול להיות ההבדל בין הימצאות מתחת לסף הקיבולת של מערכת הבריאות או קריסה שלה. כמובן, שזו רק סימולציה אחת. נדרש כעת ניתוח של המודל שיאפשר להעריך בצורה טובה יותר את השפעת אי-הציות.

מה לגבי מדיניות של בידוד האוכלוסייה בסיכון גבוה במקביל להסרת מגבלות משאר האוכלוסייה באופן שיאפשר את פתיחת מערכת החינוך והמשק באופן מיידי? כעת קבוצה 1 תהיה האוכלוסייה הכללית וקבוצה 2 תהיה האוכלוסייה בסיכון. נניח שחמישית מהאוכלוסייה בסיכון. נבחר $\beta_{11}=0.3$ ולעומתו $\beta_{22}=0.01$ שמשקף בידוד קפדני בין כל זוג אנשים באוכלוסיית הסיכון. נבחר פרמטרים שישקפו סיכויי הדבקה מסוימים בשל מגע בין האוכלוסיות, $\beta_{12}=\beta_{21}=0.05$. הגרף הבא מציג את הגודל היחסי של הנדבקים מקבוצת הסיכון בתרחיש שבו הם נמצאים בבידוד שלושה חודשים או ארבעה חודשים ולעומתו תרחיש ייחוס בו אין בידוד כלל.

Simulation of high risk/low risk gr

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

איך מתקדמים מכאן?

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

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

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

כאן זה המקום עבורכם להעיר ולשאול

בחינת נתונים עדכניים והאם סגר הוא הפתרון?

בחינת נתונים עדכניים והאם סגר הוא הפתרון?

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

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

נתונים עדכניים ומה הם משקפים?

הגרפים הבאים מציגים את מספר החולים המאומתים במספר ארצות עדכניים לתאריך 5.4, כפי שנלקחו מהאתר worldometer.com.

Data from several countries and SIR fits

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

אפשרות אחת היא שהשתנו התנאים באופן שצריך להשתקף בפרמטרים של המודל. במקרה כזה, היות והמודל מניח פרמטרים קפואים בזמן, נוצרת סטייה בין תחזית המודל לבין הנתונים. אכן, בכל המקרים הסטייה מתרחשת זמן מה לאחר שננקטו צעדים ראשונים למיגור המגפה. איני מתיימר לכמת את משמעות הצעדים שננקטו בארצות השונות באופן שיאפשר השוואה בין תחזיות המודלים. בישראל הסטייה מהאקספוננט מתרחשת שבועיים לאחר סגירת בתי הספר ובהחלט יכולה להעיד על הורדת מקדם ההדבקה בעקבות צעד זה ואחרים.
אפשרות שנייה היא שמספר החולים המאומתים אינו משקף נכונה את מספר החולים. זהו חשש סביר מכיוון שכדי לעקוב אחר העלייה האקספוננציאלית נדרשת גם עלייה אקספוננציאלית במספר הבדיקות. הסטייה מהאקספוננט החלה בישראל ב 26.3 (היום ה 24 בגרפים, כאשר יום האפס נקבע באופן שרירותי). באותו יום דיווח משרד הבריאות שבוצעו 5,300 בדיקות מתוכן $13\%$ היו חיוביות. הטבלה הבאה מפרטת את מספר הבדיקות שנדרשו כדי לעקוב אחר האקספוננט הראשוני בהנחה מקלה ש $15\%$ מהבדיקות חיוביות. לשם הבהרה, זו הדינמיקה הצפויה ללא נקיטת צעדי מנע כלשהם, ולכן השורה הזו מייצגת תרחיש מחמיר.
השורה השנייה מייצגת תרחיש מקל שיידון בהמשך.

תאריך	26.3	28.3	30.3	1.4	3.4	5.4
מספר בדיקות ליום הנדרש כדי לעקוב אחר האקספוננט הראשוני	6,965	10,896	17,046	26,667	41,719	65,266
מספר בדיקות ליום הנדרש כדי לעקוב אחר אקספוננט עם $\beta=0.19$	3,022	3,860	4,930	6,298	8,045	10,276
מספר בדיקות בפועל	6,277	5,069	5,680	8,338	7,166	7,078

האם מספר הבדיקות אכן גדל בהתאם והאם ההנחה ש $15\%$ מהבדיקות חיוביות מקלה או מחמירה בראי הנתונים?

משרד הבריאות פרסם רק ב 7.4 את נתוני הבדיקות אחרי הפסקה ארוכה מה 26.3. הטבלה מפרטת את מספר הבדיקות שנעשו בפועל, לא כולל בדיקות חוזרות לאותו אדם. ההנחה אכן מקלה – כ $12\%$ מהבדיקות חזרו חיוביות, ולא – כמות הבדיקות לא גדלה באופן שיאפשר מעקב אחר האקספוננט הראשוני.

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

איזו מהסיבות הנ״ל הובילה לסטייה של תחזית המודל מהאקספוננט?
סביר שכל התשובות נכונות במידה מסוימת. זה אחד האתגרים בניתוח מודל ובשימוש בו להבנת המציאות.

אנסה, בכל זאת, להתמודד עם השאלה. כצעד ראשון בזיהוי המקור לסטייה, אבדוק האם ניתן להסביר את הסטייה מהאקספוננט ע״י שינוי בפרמטרים של המודל בנקודת זמן כלשהי. נזכור שמעריך האקספוננט בשלבי ההתפרצות הראשונים שווה ל $\beta-1/d$ כאשר $\beta$ הוא מספר האנשים הממוצע שנשא מדביק ביום ו $d$ הוא מספר הימים הממוצע בו הוא מדבק (ימים עד להחלמה).

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

נתמקד, לדוגמא, בגרמניה. ניתן לראות שיש סטייה מהאקספוננט ביום ה 24. הנתונים לאחר יום זה עדיין גדלים באופן אקספוננציאלי אך עם מעריך קטן יותר שמתאים ל $\beta=0.19$. לאחר מספר ימים, שוב יש סטייה מהאקספוננט (קו-נקודה כחול) אך שוב הנתונים גדלים בקצב אקספוננציאלי עם מעריך שמתאים ל $\beta=0.12$. סה”כ נתוני גרמניה מתאימים היטב לפתרון של SIR בו הפרמטר $\beta$ משתנה בשתי נק’ הזמן הנ”ל. באותו אופן, ע”י שינוי מתאים בערכי $\beta$ ניתן לקבל התאמה לנתוני יתר הארצות, ראו עקומות בגרפים מעלה.

Focus on Germany

נתמקד כעת בישראל. אמנם ניתן לקבל התאמה סבירה לנתונים כאשר לאחר ה 26.3 הערך של הפרמטר $\beta$ משתנה ל $\beta=0.19$, בדיוק כמו בגרמניה, ראו עקומה כחולה בגרף המתמקד בישראל. ניתן לחזור כעת לטבלת מספרי הבדיקות ולראות אם מספר הבדיקות גדול יותר מהמספר הנדרש לעקוב אחרי אקספוננט עם מעריך המתאים ל $\beta=0.19$. התוצאה גבולית. ניתן, למעשה, לקבל התאמה טובה יותר לקצב גידול ליניארי לאחר ה 26.3 עם שיפוע שמתאים ל 530 חולים מאומתים חדשים בממוצע ביום, ראו עקומה אדומה מקווקוות. היות והמודל לא מתאר גידול ליניארי, ניתן להסביר התאמה כזו בכך שמספר החולים המאומתים נקבע ע”י מספר הבדיקות המבוצע ביום ומהווה רק חסם תחתון למספר החולים האמיתי.

Focus on Israel

האם סגר הוא הפתרון?

ניתוח הנתונים בארצות השונות מראה שעם נקיטת צעדי בידוד חברתי כאלו ואחרים, מספר האנשים שנשא ממוצע מדביק ביום פוחת באופן משמעותי ובהתאם גם קצב הגידול באוכלוסיית הנשאים. הדיון עד כה התמקד בשלב ההתפרצות ההתחלתי בו מספרי החולים מהווה שבריר מהאוכלוסייה הכוללת. נבחן עכשיו את השפעת צעדי הבידוד והסגר, כלומר את השפעת שינוי $\beta$ על מהלך המגפה כולה. תרחיש הייחוס הוא אפס התערבות, במצב כזה מגיעים לאחר כחודשיים לכ 4 מליון חולים, רובם קלים אך המיעוט שיזדקק לאשפוז יעלה בהרבה על סף הקיבולת של מערכת הבריאות. בתרחיש כזה מהלך המגפה כולה ייארך 4 חודשים עד למיגורה $98\%$ מהאוכלוסייה.

נתחיל במקרה בו מיישמים מדיניות בידוד מחמירה לתקופה מוגבלת ולאחריה חוזרים לשגרה. נבחן את ההשפעה של סגר מוחלט למשך חודשיים שממודל ע״י הפחתת הפרמטר $\beta$ ל $\beta=0$ בזמן הסגר וסגר חלקי בו $\beta=0.1$ בזמן הסגר. בשני המקרים, במהלך החודשיים של הסגר (יום 24-83) , אכן מספר החולים קטן מאוד. עם זאת, עם הסרת הסגר, אחוז המחוסנים באוכלוסייה זניח והמערכת נמצאת קרוב לאותו שיווי משקל לא יציב ממנו היא התחילה. בהתאם המגפה מתפרצת שוב (ראו עקומה אדומה) עם מאפיינים כמעט זהים למהלך המגפה בתרחיש הייחוס. המספר הכולל של החולים לאורך מהלך המגפה נשאר כמעט ללא שינוי. הסגר החלקי מיטיב מעט את המצב בהיבט של מספר החולים בשיא ומספר החולים הכולל. זה נובע מכך שמצטברים מחלימים בעלי נוגדנים שיכולים לתרום לחיסוניות העדר, אך השיפור מזערי. בראייה כוללת, השפעת סגר, אפילו מוחלט, היא פשוט דחייה של מהלך המגפה.

נבחן כעת תרחיש אחר בו הפרמטר $\beta$ קטן לאורך תקופה ארוכה, אך לא בשיעור חד כמו בסגר חלקי. נקיטת צעדים שמפחיתים לאורך זמן מ $\beta=0.3$ ל $\beta=0.19$ מביאים לשיא של כ $2.6$ מליון חולים לאחר כשלושה חודשים וחצי מההתפרצות ההתחלתית. עד למיגור המגפה $91\%$ מהאוכלוסייה תחלה. זו ירידה במספר החולים הכולל שיחלה. מכך נגזרת גם ירידה בסך מקרי המוות כתוצאה מהמגפה. עם כל זאת, מספר החולים בשיא גבוה ונגזרים ממנו צורכי אשפוז שהם עדיין הרבה מעבר לסף הקיבולת של מערכת הבריאות. הקטנת $\beta$ לאורך זמן במקביל לשחרור המשק היא אתגר, שיחייב אימוץ נורמות התנהגות חדשות ע”י האוכלוסייה הרחבה.

השפעה של שינוי בטא

איך אפשר לרתום את מנגנון חיסוניות העדר?

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

בפוסט הבא, אציג הרחבות של מודל SIR שמאפשרות לבחון השפעה של בידוד סלקטיבי ואשתמש בהן כדי לענות על מספר שאלות, כגון: עד כמה מהלך כזה יכול להיות אפקטיבי? כיצד ליישם באופן אופטימלי מדיניות של בידוד סלקטיבי? עד כמה מדיניות כזו תלויה באחוז ציות להנחיות – מה קורה אם רק $70\%$ או אפילו $95\%$ מהאוכלוסייה מקפידה על ההנחיות?

כאן זה המקום עבורכם להעיר ולשאול

מהלך המגפה

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

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

מהלך המגפה הכולל כפי שמתאר המודל

נסתכל על המשוואה ל $I(t)$ הרשומה בעזרת מקדם ההדבקה ונבחן אותה במקרה בו $R_{0}>1$,

$$ I^{\prime}(t)=\beta\left(S(t)-\frac{1}{R_{0}}\right)I(t) $$ בהתחלה $S(t)\approx1$ומכיוון ש $R_{0}>1$ נובע ש $I^{\prime}(t)>0$. העלייה במספר הנדבקים נמשכת כל עוד $S(t)-\frac{1}{R_{0}}>0$. תהליך ההידבקות הוא תהליך של מעבר אוכלוסייה מקבוצת המועדים לקבוצת הנשאים. ככל שעולה מספר הנשאים, כך יורד מספר המועדים. בנקודת זמן $t_{m}$ כלשהי הגודל היחסי של קבוצת המועדים מגיע לסף $S(t_{m})=\frac{1}{R_{0}}$. זוהי נקודת המפנה בה העלייה במספר הנשאים נעצרת ($I^{\prime}(t)=0$) ומגיע לשיא. מרגע זה, מספר הנשאים מתחיל לרדת בהדרגה עד שהוא מגיע לאפס. בשלב זה המגפה מסתיימת, והמערכת מגיעה לשיווי משקל חדש. זהו המהלך הטבעי של המגפה המתוארת ע״י המודל ללא נקיטת צעדים כלשהם להדברתה. מהלך המגפה מומחש בגרפים של משתני המערכת (מבטיח השוואה לנתונים בפוסט הבא). בפרט, הגרף של מספר הנשאים הוא מה שמתייחסים אליו לפעמים כאל ה׳עקומה׳.

Conceptual Dynamics of SIR

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

מספר החולים בשיא המגפה

אתחיל משאלת חישוב מספר הנשאים בשיא המגפה – קיבלנו שבשיא המגפה מספר המועדים הוא $S(t_{m})=\frac{1}{R_{0}}$. כעת אנסה למצוא קשר בין גודל קבוצת המועדים לגודל קבוצת הנשאים. זהירות, קטע מעט טכני לפנייך – אם נחלק את המשוואה עבור $I(t)$ במשוואה עבור $S(t)$ נקבל $$ \frac{I^{\prime}(t)}{S^{\prime}(t)}=\frac{\beta S(t)I(t)-I(t)/d}{-\beta S(t)I(t)} $$

נוכל להפריד אגפים ולקבל $$ I^{\prime}(t)=\left(\frac{1}{\beta dS(t)}-1\right)S^{\prime}(t) $$

כעת ע״י אינטגרציה של כל אגף ושימוש בטרמינולוגיה של מקדם ההדבקה נקבל את הקשר הרצוי בין ל $I(t)$ ל $S(t)$ $$ I(t)=1-S+\frac{1}{R_{0}}\log S $$

בפרט, בשיא המגפה, אחוז הנשאים באוכלוסייה הוא

$$ I(t_m)=1-S(t_m)+\frac{1}{R_{0}}\log S(t_m)=1-\frac{1}{R_{0}}-{\frac{\log R_{0}}{R_{0}}} $$

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

אקדיש רגע להתייחס לנוסחא המתקבלת. מספר הנשאים בשיא תלוי רק במקדם ההדבקה $R_{0}$. המשמעות היא שעל פי המודל הדרך היחידה לשטח את העקומה היא להפחית את מקדם ההדבקה. זו אחת הסיבות שמקדם ההדבקה תופס כל כך הרבה תשומת לב לאחרונה. המודל ׳רואה׳ רק את מקדם ההדבקה ואדיש לדרך בה פועלים להפחיתו. לדוגמא, צמצום אינטרקציות בתוך האוכלוסייה או איתור יעיל של נשאים חדשים ובידודם. מהו ערך המטרה שיאפשר למערכת הבריאות להכיל את המגפה? תלוי במערכת הבריאות. הגרף הבא מציג את אחוז הנשאים באוכלוסייה כפונקציה של $R_{0}$כאשר האזור האפור ממציין את טווח הערכים של מקדמי ההדבקה בארצות שונות כפי שחושבו בפוסט הקודם. בטווח מקדמי ההדבקה שנמדדו בשלבים הראשונים של ההתפרצות, המודל מנבא שבשיא המגפה כ $ 45\%-40\%$ מהאוכלוסייה תהיה חולה.

Graph of I(tm) as a function of R0

מהו המנגנון שבולם את התפרצות המגפה?

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

אפשר לשאול גם כמה יחלו עד לשיא המגפה, כלומר לכלול גם את מספר המחלימים. במקרה זה התשובה כבר זמינה, $$ I(t_m)+R(t_m)=1-S(t_m)=1-\frac{1}{R_{0}}$$ בטווח מקדמי ההדבקה שנמדדו מגיעים ל $ 79\%-74\%$ מהאוכלוסייה. זה האחוז הנדרש כדי שחיסוניות העדר תבטיח דיכוי טבעי של המגפה. באפידימיולוגיה המספר הזה נקרא סף חיסוניות העדר. ייתכן שזה המקור להערכות ש $ 70\%$ מהאוכולוסייה יחלו בקורונה, אם כי המספרים האלו מתייחסים רק לאחוז האוכלוסייה שנדבק עד להגעה לשיא המגפה ולא לשלב בו המגפה דועכת בהדרגה. בטווח הפרמטרים הנ״ל, בתקופה הדעיכה ההדרגתית של המגפה עוד כרבע מהאוכלוסייה תחלה עד שפתרון המודל יתכנס לשיווי המשקל בו רק $ 2\%-1\%$ נשארים מועדים.

אחת המסקנות העיקריות מניתוח המודל עד כה היא שהדרך היחידה להבטיח את בלימת ההתפשטות של הנגיף לטווח ארוך היא לעבור את סף חיסוניות העדר. אפשר לשאול האם יהיה אפקטיבי להוריד את סף חיסוניות העדר ע״י הורדת מקדם ההדבקה ארוך הטווח. לדוגמא ע״י הטמעת נורמות התנהגות חדשות באוכלוסייה. זו למעשה שאלה שעוסקת ברגישות המערכת לשינויים – האם מספיק, לדוגמא, ש $ 80\%$ מהאוכלוסייה יקפידו על נורמות התנהגות המותאמות לעידן ה(פוסט-)קורונה כדי להפחית את סף חיסוניות העדר? אחד הפוסטים הבאים יעסוק בצעדי מנע קצרי טווח וארוכי טווח.

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

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

כאן זה המקום שלכם להעיר ולשאול

שלב ההתפרצות הראשון

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

בחלק הראשון הצגתי את פיתוח מודל SIR – מודל המחלק את האוכלוסייה לשלוש קבוצות: קבוצת המועדים (Susceptibles), הנשאים (Infected) והמחוסנים ,(Recovered) ומגדיר את האופן והקצב שבו אנשים עוברים מקבוצה לקבוצה, לדוגמא ע״י הידבקות או החלמה. לאחר מכן הצגתי ניתוח ראשוני של המודל והראיתי שמגפה יכולה להתפשט רק כאשר מקדם ההדבקה או מספר האנשים שנשא מדביק בממוצע טרם החלמתו גדול מאחד. בפרמטרים של המודל, התנאי מתרגם ל $R_{0}=\beta d>1$ . הניתוח הנ״ל התמקד במצבי שיווי המשקל של המערכת ובחן את היציבות שלהם, כלומר את התנהגות המערכת כאשר היא מוסטת קלות משיווי המשקל. כעת אתמקד בשלב הראשון של המחלה – כאשר המערכת מוסטת משיווי המשקל והמגפה מתחילה להתפרץ.

שלב ההתפרצות הראשון

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

$$S(t)=1-\frac{10}{9,000,000}\approx0.999998 $$

גם כאשר יש 10,000 חולים, הגודל היחסי של קבוצת המועדים הוא קרוב מאוד לאחד, $S(t)\approx0.998$. לכן, בשלב הראשון וכל עוד מספר החולים קטן מאוד ביחס לגודל האוכלוסייה, אפשר להשתמש בקירוב $S(t)\approx1$ כדי לחקור את המערכת. המשוואה המתקבלת ל $I(t)$ היא $$ I^{\prime}(t)=(\beta-1/d)I(t) $$ בזכות הקירוב קיבלנו משוואה פשוטה יותר שניתן לפתור באופן מפורש. הפתרון של המשוואה הוא אקספוננט $ I(t)=I(0)e^{(\beta-1/d)t} $. המודל מתאר אם כן גידול אקספוננציאלי בשלבים הראשונים של המחלה.

אתייחס כעת לנתונים של כמות החולים המאומתים במספר ארצות. לשם נוחות הגרפים מציגים את הגודל המספרי של קבוצת הנשאים, כלומר את $N\cdot I(t)$ כאשר $N$ מייצג את גודל האוכלוסייה. האיקסים הכחולים הם נתוני כמות חולים אשר נלקחו מהאתר www.worldmeters.info, ואילו העקומה האדומה בכל גרף היא האקספוננט המתאים ביותר לנתונים בשבועות הראשונים.

Graphs of number of infected people in Corona as a function of time in various countries

בכל ארבעת הגרפים השלב ההתחלתי אכן מתואר היטב ע״י גידול אקספוננציאלי. הנתונים העולמיים מצביעים על זמן החלמה ממוצע של שבועיים ובהתאם אבחר $d=14$. כעת ניתן לחלץ את $\beta$ ממעריך האקספוננט המתאים ביותר לנתונים. הערכים המתקבלים נעים בטווח של 0.28-0.32 בארצות האלו, כאשר יש ארצות עם $\beta$ נמוך יותר (איטליה, לדוגמא, עם $\beta\approx0.25$). בתקשורת נוטים להציג את ערך המעריך בעקיפין דרך מספר הימים בהם מספר החולים מכפיל את עצמו. במקרים המוצגים קבוצת החולים מכפילה את גודלה כל 2.82-3.27 ימים.

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

בחלק הבא אתמקד במהלך המגפה המלא כפי שמתאר המודל ואנסה לשפוך אור על חלק מהשאלות.

נ.ב.

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

כאן זה המקום שלכם להגיב ולשאול

זמן	0	\( \tau\)	\( 2\tau\)	…	\( 6\tau\)	\( 10\tau\)	\( 15\tau\)
נדבקים בזן המוכר \( I_1(t)\)	1000	1200	1444		3000	6200	15,400
נדבקים בווריאנט הבריטי \( I_2(t)\)	10	18	32		340	3600	67,500
אחוז נדבקים בווריאנט הבריטי \( \frac{I_2(t)}{I_1(t)+I_2(t)}\)	0.9%	1.5%	2.2%		10%	36%	81%