2020年8月25日火曜日

「ヨシッ、解った !?」検査の確率 - 夏休み「明け」の自由研究

残暑厳しい折、お見舞い申し上げます。

自由研究のネタ探しに来られたお子さん、今年の夏休みはコロナウィルス感染症蔓延のため、早々にお休みとなって、夏休み本番はかなり短いものとなってしまった様です。

有名どころお医者さんを中心に、コロナウィルス感染のPCR検査は不要という論調でありますが、その辺りを深く考えるため、PCR検査などの検査について、少し調べてみました。

PCR検査などには、「感度」、「特異度」というパラメタがあるらしい。

+ 「感度」  実際に感染している人の検出の率 (a としておきます)
+ 「特異度」感染していない人の検出の率 (b としておきます)

感度は、実際に感染している人が10人いたとして、その内、どれだけの人が「感染している」事を検出できるか、という割合なんだそうです。
例えば、感度が70%とすると、10人の感染者の内、70% = 7人しか検査で検出できない、ということです。
これは、10人の感染者がいた場合、感度 70% では、3人の取りこぼしがある勘定です。この取りこぼしの事を「偽陰性」と言い、偽陰性の割合は「偽陰性率」ですから、

  感度 + 偽陰性率 = 100 %

になります。

特異度は、感染していない人の内、どれだけの人が「感染していない」事を検出できるか、という割合だそうです。
特異度が90%とすると、10人の感染していない人の内、検査によって感染していないと検出できる人は、90% = 9人、となります。
特異度90%ならば、10人の内、1人が「取りこぼし」となり、その取りこぼしの人は、実際には感染していない(陰性)にも関わらず、検査では陽性=「偽陽性」となってしまいます。

  特異度 + 偽陽性率 = 100 %

感度、特異度、などと書いていると、キーを打つ回数が増えてしまうので、ここからは感度の代わりに a 、特異度の代わりに b と書きます。また、パーセント表記は、これも混乱しやすいので、100 % を 1 とした表記で進めていきます。

  感度  (a) + 偽陰性率 = 1
  特異度 (b) + 偽陽性率 = 1

  偽陰性率 = 1 - 感度  (a) = 1 - a
  偽陽性率 = 1 - 特異度 (b) = 1 - b

となります。

ここまでが、検査についての予備知識です。
今般、流行しているコロナウィルス感染症ですが、感染の割合が、よく判っておりません。
そこで、取り敢えず、感染している人の割合を、0.1 % としておきます。

この割合で言えば、例えば、10000人ならば、感染者は10人くらいになります。感染していない人は、10000-10 = 9990 人です。

この10000人に対して、感度 a = 0.7, 特異度 b = 0.9 の検査を行うと、

感染者 10 人 :
  検査陽性 = 10 * 0.7 = 7 人、  検査陰性 = 10 * 0.3 = 3 人

非感染者 9990 人 :
  検査陽性 = 9990 * 0.10 = 999 人、 検査陰性 = 9990 * 0.90 = 8991 人

となるのだそうです。

実際に感染していない人でも、検査で「陽性」となってしまう人が999人も出てしまう、PCR検査はむやみに行うべきではない、というのが、有名なお医者さんの意見だそうです。

有名どころのお医者さんの説明は、つぎのサイト様の記事で、色々と紹介されております。

cf. ベイズの定理を悪用し、コロナウイルスPCR検査の有用性を否定する医師達 - 臨床獣医師の立場から
https://tatsuharug.com/abuse-bayes?fbclid=IwAR0WTAjomzGlZWvilmT_XalTqlNNlJ3Y4R-nPJKu7-LFvl9sdXEFNUiRcn4#i-8

上記の記事によりますと、有名どころのお医者さんの多くは「PCR検査の特異度は、90%、99% にとどまる」としているのだそうですが、
記事を書かれている方は、「PCR検査の特異度は、(90%、99% ではなく) 99.99 % 以上」と述べておられます。

素人目には、99%も99.99%も、違いがある様には思えないのですが、果たして、実態はどうなのか ?

非感染者で検査陽性となる人については、つぎの計算式で計算できます。

  (非感染、検査陽性数) = N*(1-p)*(1-b)

  N ; 検査総数
  p ; 感染者の存在確率
  b ; 検査「特異度」

感染者の存在確率、特異度が低いと、 この値は大きくなります。
そして、特異度の値の代わりに「偽陽性率」(= 1-b)でみると、より判りやすい。

+ 特異度 b = 0.90 (90 %) ならば、偽陽性率 1-b = 0.10
+ 特異度 b = 0.99 (99 %) ならば、偽陽性率 1-b = 0.01
+ 特異度 b = 0.999 (99.9 %) ならば、偽陽性率 1-b = 0.001

特異度の精度が向上するだけで、偽陽性率が劇的に下がります。
偽陽性率が劇的に下がれば、偽陽性者の数は十分下がり、有名どころのお医者さんが言っている「検査大パニック」は起こらない。
上記のサイト様では、この事を丁寧に述べております。

PCR検査、今は検査数が大分増えてきて、感染者も多く確認出来ております。
しかしながら、PCR検査、多くの場合「自費」で行うらしく、健康保険も使えないので、かなりの額を自腹でやらんとアカンとか。
また、市中では相変わらず「感染リスク」があるので、検査で陰性と結果が出ても検査実施機関では「陰性証明書は出しません」。
海外では、PCR検査が安価に実施出来るというのに、本邦では、健康保険の適用外で「自費」。しかも、有名どころのお医者さんは「タダの風邪なんだから、家で寝てれば治るっぺ。だから、検査もせんでエエ。下手に検査で来られても、検査大パニックだかんネ」というつもりらしい。

こうした(特異度の余り高くない)検査について、面白い話題があります。

cf. 【受験数学】条件付き確率の公式とイメージを徹底解説!!【確率】(例題つき) - hmorinari's diary
https://hmorinari.hatenablog.com/entry/2019/01/18/224939

このサイト様によると、1回の検査で陽性の結果が出ても、感染している確率は低いものだが、その検査陽性者集団に再検査すると、さらに確度が高くなると見込まれるらしく、「感染している確率が低いとは言え、検査自体が無意味ではありません」と述べています。

感染者の存在確率 (「事前確率」というそうです)を p としておきます。
サイト様の説明によると、「検査を受け、陽性という結果になった人が、ホントに感染しているのか」その確率 p' は、

  p' = (p*a) / ( (p*a) +  (1-p)*(1-b) )

で計算できるとのこと。
但し、
  p' ; 検査の結果、陽性だった人がホントに感染している確率
  p ; 感染症に罹患している、平均確率
  a ; 検査「感度」
  b ; 検査「特異度」
です。

サイト様の例では、
  p = 0.001 (0.1 %), a = 0.95 (95 %), b = 0.98 (98 %)
で計算しているので、検査を受け、陽性の結果になった人(サイト様では「A君」としております)が、ホントに感染している確率は

  p' = 0.04538 ... (4.538 %)

ですが、これは、同時に、検査で陽性の結果になった人々の感染確率=(検査陽性集団の)事前確率、とも言えるものです。
そこで、A君を含めた検査陽性者集団を再検査すると、

  p'' = 0.693... (69.3 %)

になります。再検査で陽性になった人々に、再再検査をすれば、

  p''' = 0.9907... (99.07 %)

という確度になります。

こうした、ちょっとした計算式の繰り返し適用は、フツーの電卓で行うには少々手間。一方、表計算ソフトで行うのは「牛刀を以って...」という例えの様に大掛かりではあります。
そこで、高機能電卓を用いるのは、どうでしょう。

感度 a と、特異度 b について、ある関係 (a+b > 1)を満たしていれば、多重検査による確度の加速は達成できます。

今度は、つぎの例でやってみました。
  p = 0.001 (0.1 %), a = 0.7 (70 %), b = 0.90 (90 %)

   0 回め   0.0069582
   1 回め   0.0467557
   2 回め   0.2555886
   3 回め   0.7061764
   4 回め   0.9438953

ここまで特異度が低い (偽陽性度が高い)場合、確度の加速を目指すとなると、再検査過程が多数回必要になります。実際、医療の現場で、こんな悠長な事をしているのだろうか ? やはり、多くの有名どころのお医者さんが申しておる、特異度 90〜99 % ちう値は、実用的ではないのだろうと。
実際、検体を放り込むと、あとは自動で検査をしてくれる医療機器とかあるそうで、こうした機器に採取した検体を放り込んでしまえば、検査過程でのコンタミ (検査過程での検体の汚染、コンタミネーション)も無さそうです。コンタミがなければ、特異度 99.99%以上、とか達成できそうです。

さて、実際はどうなのでしょうか ?

2020年8月4日火曜日

(ヤメては居りませんのヨ)

何か、書くことも思いつかず、まずは消息ポストです、ハイ。
申し訳ない。

以上

2020年4月12日日曜日

新型肺炎禍 ...

新型肺炎の蔓延を受け、なるべく外出をしないように、との事です。

海外では、こんな話題があります。

Texas Instruments、COVID-19で自宅学習を余儀なくされている学生を支援するため、グラフ電卓アプリ「TI-Nspire CAS for iPad」などを一時的に無料に。 - APPL CH
https://applech2.com/archives/20200321-texas-instruments-ti-nspire-cas-for-ipad-mac-free.html

TIはアメリカでグラフ電卓のシェア#1なので、こうした事をしております。
CASIOも負けずに、こうした取り組みをしております。

基本操作から例題による活用まで! 動画でわかるグラフ関数電卓 - casio.jp
https://web.casio.jp/dentaku/fxcg50/movie/?topics

こうした企業努力は素晴らしいものです。

さて、当初は「検査など不要」としていた厚労大臣でしたが、今日では、多くの自治体で新型肺炎の検査をしているとか。

東京都のバヤイ、2020/04/12の感染者数は197人となっております。

cf. 都内の最新感染動向
https://stopcovid19.metro.tokyo.lg.jp/

しかし、検査した人数は344人、との事。エエッー、少ないじゃん !!

当方、高機能電卓を使っておるので、簡単な計算をしてみましたヨ、ウヒョヒョ。

検査数に対する感染者の割合、という事で、こうした場合には「母比率の推計」という機能を使います。検査総数 n=344, 感染者数 x=197, 信頼区間 C=0.99 として計算すると、

hat p = 0.57267
0.5039 < p < 0.6413

となります。自覚症状があって検査に来ている人の、かなりの割合が陽性となっているらしい。

東京都の検査数、実際にはかなり少ないものですが、それもそのはずで、

(注)医療機関が保険適用で行った検査は含まれていない

なのだと。一体、誰を調べているのか ? 自費で検査した人だけなのか ?

医療機関の検査実数が上がってこない状況で、医療現場の混乱は伝わって来ない。
現実には、感染爆発が進行しているのかも知れません。

2020年3月14日土曜日

今日はホワイトデーなので (?)、久しぶりにHP Primeのコード

まずは、こちらのtweetを御参照

cf,
https://twitter.com/RR_Inyo/status/1217489174128758784

z' = cos(Z) の繰り返しで、こんなグラフィクスが描けるそうで、HP PRIMEでもやってみましたヨ。
まだ、実機を持っていないので、今回もSIMで作業ですが、以下のコードを御利用下さい。

// title : iter_cos.ppl - iterative complex cosine graph for HP Prime
// begin : 2020-02-16 20:37:33 

//
pixVal(z, c)
BEGIN
  LOCAL  cnt := 0;

  //  この判断式がミソです
  //  Mandelobrot, julia setでは、abs(z)<4 .0="" br="">  //  このグラフィクスのバヤイはこんな風でした
  WHILE (ABS(z) <= 50) AND (cnt < 20) DO
    // modify this form to enjoy
    z := COS(z);
    cnt := cnt+1;
  END;

  //  color value : 0 for black, 16777215 for white
  RETURN RGB(1*cnt/20*255, 0.75*cnt/20*255, 0.25*cnt/20*255);

END;

EXPORT itercos()
BEGIN
  // Clean the screen (G0):RECT();

  LOCAL dx, dy, z, xp, yp; 
  LOCAL xmin, xmax, ymin, ymax;
  LOCAL pixStat;

  //  screen range
  xmin := -10;
  xmax :=  10;
  ymin :=  -5;
  ymax :=   5;

  dx := (xmax-xmin)/320;
  dy := (ymax-ymin)/240;

  // we loop over each pixel
  // Of the HP Ptime screen 
  FOR yp FROM 0 TO 239 DO
    FOR xp FROM 0 TO 319 DO
      z := (xmin+xp*dx, ymax-yp*dy);
      pixStat := pixVal(z, Z0);
      PIXON_P(xp, yp, pixStat);

    END;
  END;

  //  key wait loop
  REPEAT UNTIL GETKEY() == -1;
  FREEZE;

END;



Numworks, mandelbrotのコードを改変したとの事で、単純に式を置き換えたのですが、それだけでは、こうしたグラフィクスにはなりませんでした。
色々といじくり回し、最終的に繰り返しの判断部分をいじる事で、Tweetの様なグラフィクスを得た次第。

最近、HP Prime G2の価格が大分下落してきております。

cf.
https://www.amazon.co.jp/gp/offer-listing/B07HF6RXGG/ref=dp_olp_new?ie=UTF8&condition=new

G1との価格差も狭まってきました。Sentaro様曰く、G2は動作軽快でメモリ増量。ウーム。

一方、G2は、新しいファームがあるそうです。

cf. HP Prime Beta Software
https://www.hpcalc.org/prime/beta/

ハードウェアの違いがあるので、ファームも更新されたらしい。

そういや、こうしたコードを挙げておきながら、HP Prime実機への転送方式、気にした事がありませんでした。どなたか、ご教示ただけると幸いです。やはり、HP Conn-Kitしかないのかなァ ?

2020年1月2日木曜日

С Новым годом 2020 !


昨年もお世話になりました。今年も宜しくお願い申し上げます。

今年はロシア語で「謹賀新年」です。web時代で、簡単に見つかりますネ。

昨年末の慌ただしい中、性懲りもなくこんな計算をしておりました。

cf. ぽてん (@potenten87) 様 tweet
数学ができないとWi-Fiに繋げないホテルがロンドンにあるって。
https://twitter.com/potenten87/status/1190161683638960130



〽ロンドン、ロンドーン ! 楽しいロンドン愉快なロンドン、ロンドンロンドーン !!
(歳がバレますなァ)

この式、計算すると円周率piになるそうで、電卓でもって早速やってみた次第。



Numworks web simlatorの結果



HP Prime β版simulator




数値積分をやってみると、πに近い値になるようです。

Numworksのweb sim (https://www.numworks.com/simulator/) はJavascriptで実現されているそうで、もしかすると実機とは異なる答を出すかも知れません。

また、HP Primeのシミュレータはβ版なので、結果はこんな値でした。
XCASが使えるのですが、やり方、よく判らんですよ。申し訳ない。

量販店に行き、店頭のfx-CG50でも試してきました。
すると、CG50ではズバリπを返しますネ ! これはオモシロイ。

ただ、角度モードがRADになっている場合に限るようです。

それは、数値積分を愚直にやっていて、RADモードにした際の計算誤差がπに近いからなのだろうか、と想像します。CASIOの数値表記機能がここで威力を発揮したのでしょう。

計算の結果がπになる、ちうのは、上記のTweetをみていくと、色々出ております。奇関数、偶関数んて、すっかり忘れておりましたヨ。オツムの足りない当方なので、久しぶりに筆算をして、ようやく納得がいった次第です。今時の若い人も、こんな計算をして大学受験に備えておるのですが、大学出てしまうと、すっかり忘れてしまう人も多いかも ? (いや、ソレは当方だけなのか ?)