電卓で「濃厚接触」 - dentagion

HP Primeで、簡易な「感染シミュレーション」をやってみました。

#  あらまし

簡単な感染シミュレーションです。

100人程度のごく小さな集団にて、感染者が、どの様に増えていくのか、それを数値実験するものであります。
リスト L1 に100人の状態を設定しております。

リストの要素(セル)1つが1人の「感染状態」を保持していて、
  • 0 - 未感染
  •  1〜3 - 感染
  •  4 - 発症/隔離
としています。

リストのセルを一通り調べ、
  1.  セルの状態が4ならば「発症/隔離」として、特に何もせず
  2. セルの状態が1〜3ならば「感染」として、感染の状態数を1増やして、更新
  3.  セルの状態が0ならば「未感染」で、両隣のセルを調査
     左右のセルが感染していたら、乱数によって「rate」を上回った場合に、当該未感染セルを「感染」として、状態数1 にする

という処理を実施します。

最初に、リストの中央に感染セル(状態数 1)を設置、あとは、次の流れでシミュレーションを進行します。
  1. セルのランダム移動
  2. 感染処理
  3.  未感染(Suscept), 感染(Infection), 発症/回復(Recover) の集計
この処理を30回ほど実施、経過を C1 (インデクス), C2 (未感染数), C3 (感染数), C4 (回復数) の時系列リストとして保存して終了します。

#  使い方

プログラム「dentagi.ppl」を読み込んで、実行して下さい。
最初に「rate」と、数値入力のテキストボックスが表示されるので、0から1の範囲で、適当な数字を入力してやります。すると、リスト C1, C2, C3, C4 に数値を保存し、実行を終了します。

HP Primeでは「2変数統計」機能で、リストのscatter plotを実施できますから、こうして出来上がったC1, ... C4 のリストデータを用いて、感染者数の推移をグラフに出来る、という具合。
  • 未感染のグラフを描くには、C1 (X軸), C2 (Y軸)
  • 感染のグラフを描くには、 C1 (X軸), C3 (Y軸)
  •  発症/回復のグラフを描くには、C1 (X軸), C4 (Y軸)
と設定しておきます。処理ループは30回なので、X軸の範囲は、0〜30 の程度に設定しておけばよさそうです。

図は、rate=0.5の例で、100人程度の「濃厚接触空間」でも、SIRモデルっぽい曲線が描かれます。



一度、プログラムを読み込んで置けば、後はHOMEボタンを押しdentagi(0.5) の様に打ち込んでやれば、様々な感染レートでのシミュレーションが実施できます。

感染レートが低い場合には、未感染の数がそれなりに増える場合がある、など、様々な事が観察されます。

色々とお試し戴き度。

以下、リストだよん。

// title : dentagi.ppl - micro contagion simulation for HP Prime
// begin : 2020-04-27 00:56:02 
// note  : 

//  L0 ; cell space
//  L1 ; cell space next stage
//  C2 ; Suscept, C3 ; Infect, C4 ; Recovre

//  rate ; infect Rate [0, 1)
EXPORT dentagi(rate)
BEGIN

  LOCAL  c, i, j, k, len, p1, p2, tmp;
  LOCAL  idx, sus, inf, rec; 

  //  cell space length
  len := 100;

  //  cell array initialize
  L1 := MAKELIST(0, X, 1, len, 1);
  L2 := MAKELIST(0, X, 1, len, 1);
  C1 := {};
  C2 := {};
  C3 := {};
  C4 := {};

  //  place first infected cell
  L1(len/2) := 1;

  // contagion game loop
  FOR k FROM 1 TO 30 DO

    //  random walking
    FOR i FROM 1 TO 30 DO
      p1 := RANDINT(1, len);
      p2 := RANDINT(1, len);
      tmp := L1(p1);
      L1(p1) := L1(p2);
      L1(p2) := tmp;
    END;

    //  infection part
    L2 := L1;

    FOR i FROM 2 TO len-1 DO
      //  suscept cell is infected from right infected cell
      IF L1(i) == 0 AND 0 < L1(i+1) < 4 AND RANDOM() < rate THEN
          L2(i) := 1;
      END;
      //  suscept cell is infected from left infected cell
      IF L1(i) == 0 AND 0 < L1(i-1) < 4 AND RANDOM() < rate THEN
          L2(i) := 1;
      END;
      //  infected cell stage number increase
      IF 0 < L1(i) < 4 THEN L2(i) := L2(i)+1;  END;
    END;

    L1 := L2;

    //  collect S,I,R 
    sus:=0;
    inf:=0;
    rec:=0;

    FOR i FROM 1 TO len DO
      c := L1(i);
      IF c == 0 THEN 
        sus := sus+1;
      END;
      IF 0 < c < 4 THEN
        inf := inf+1;
      END;
      IF c == 4 THEN
        rec := rec+1;
      END;
    END;

    //  add SIR summa to C2, C3, C4
    C1 := CONCAT(C1, k);
    C2 := CONCAT(C2, sus);
    C3 := CONCAT(C3, inf);
    C4 := CONCAT(C4, rec);

  END;

END;

コメント

このブログの人気の投稿

マニング公式の計算

Confidence interval for a proportion を計算するゾ !