HPO機密日誌

自己をならふといふは、自己をわするるなり。

対話型データベースとしてのG Suite、スプレッドシート

ある方からスプレッドシート特有の関数、機能について手ほどきを受けた。importhtml、queryとarrayformuraだ。importhtmlでホームページのtable要素、list要素をリアルタイムで取り込むことができる。queryは取り込んだ表をフィルターや、ソート、カテゴリー化することができる。arrayformuraは、これまで何百行もコピペしていたvlookupや、text編集の関数を一行で必要なすべてのセルに適用させることができる。いままでの関数と大きく違うのは、これからの関数は書き込んだセルだけでなく、データに応じた多数のセル上に参照可能な形で展開できることだ。

結果はものすごい(私にとっては)。

実際に使うに当たって、かなりこのサイトを参考にさせていただいた。

tonari-it.com

tonari-it.com

tonari-it.com

こんなグラフの下材料となるデータの加工もかなり簡単にできる。Pivotテーブルがあるじゃないかという人もいるだろうが、どうもデータの取り込みだったり、数があわなかったり、私にとっては使いづらかった。

「コツ」もいくつかあるので、あとで多少補足する。

Covid19大流行の国とそうでもない国(個人的メモ)

@J_Sato さんのBCG仮説が次第に注目されてきている。私も力任せに、人口三千万人以上の国のCovid19感染状況とBCG接種状況を集計してみた。

docs.google.com

ソースは、Worldmeterと教えていただいたこちら。

結果、BCG接種国は圧倒的に人口百万人当たりのCovid19関連死者数が少ないことを改めて確認させていただいた。

f:id:hihi01:20200624211903p:plain

ちなみに、このグラフの各カテゴリー毎の国の数、そして人口百万人当たりのCovid19死亡者数。

カテゴリー 国の数 平均百万人当たりCovid19死者数
BCG接種国 31 13人
例外(後でBCG接種国と確認 5 12人
南米のBCG接種国 3 212人
現在BCG非接種国 8 384人

地図で見るとこんな感じ。

f:id:hihi01:20200624212644p:plain

http://www.bcgatlas.org/

グラフを見て、国の数を考慮に入れるとBCG接種国とそれ以外で全くCovid19の感染状況が違うように見える。もっとも、人口三千万人以上の国だけを抽出しているし、南米のブラジル、ペルー、メキシコの大流行を説明出来ない。逆に言えば、大流行かつ大規模の死者数を出しているのは北米と欧州に限られるようにも見える。

なにがこれらの国々の状況を分けたのか、ぜひ医学的な根拠を明らかにしていただきたいと。


■追記

教えていただいた世界銀行等のデータを元にG Suiteのスプレッドシートのマップグラフに表現してみた。人口が少ないと特異な部分があるように思えるので人口一千万人以上の国をフィルターした。円グラフの大きさが各国の人口百万人当たりの死亡者数、色がBCGの推測接種率。Worldmeterと世銀データで国名の表現がS. KoreaとRepublic of Koreaとか、Democratic Republic CongoとD.P.Cなど整合していないのに苦労している。

f:id:hihi01:20200627062620p:plain

https://docs.google.com/spreadsheets/d/1y0JDTPBkD_MV3yRu3h3nO2RZp68rGB4NlvXupDx3piA/edit?usp=sharing

こうみると、BCGとの関連性もさることながら百万人当たり死亡者数が100人を越える国の方が少数派なのではないかと見える。日本を含む「軽症」国が特殊なのではなく、北米、欧州の方が特殊なのではないだろうか?さらには、各国で対策の強弱は相当に異なるはずなのに、BCG接種国の百万人当たりの死亡者数の分散は小さいように思える。南米のブラジル、ペルー当たりの例外は確かにある。偉い方にぜひ解明していただきたい。

BCG要素を外して純粋に人口百万人当たりの死亡者数で色分けしてみた。

f:id:hihi01:20200627065250p:plain


■追記の追記 20200630

@J_Satoさんにきちんと死亡者データを取っているサイトを教えていただいた。色の濃淡が割と分かりづらいが、大流行の国とそれほどでもない国が色分けされている。

f:id:hihi01:20200630174815p:plain

https://ourworldindata.org/grapher/total-covid-deaths-per-million?year=2020-06-29

実は、リンクを教えていただいたからかなりの時間、地図データにたどり着けなかった。「タブ」に「map」とあるのでそこをクリックするとでてくる。時系列のアニメーションにもなるので興味深い。

「超勉強力」

やはり、中野信子先生はすばらしい。好きなことを好きだと言える、好きなことをここまで極められるお姿がすばらしい。

で、とうとう買ってしまった。

「超」勉強力

「超」勉強力

中野先生は、以前から崇拝していた。

hpo.hatenablog.com

好きなことを極めることが一番の勉強法だと断言されているのことにもうツーサムスアップ!

全然レベルは違うが、まわりから幼いころこら「あんたは人の気持ちがわからない」、「もっとばからにならないと」とか私ですら言われてきた。中野先生のレベルでは更にいろいろなコンプレックスを抱かせるようなまわりの言動もあったのだろうと端々からうかがえる。

厚労省大規模抗体検査結果の発表

待ちに待った発表があった。結果はかなり予想外ではあったが。

厚生労働省は、今月1日から7日にかけて人口が一定規模ある地域のうち、10万人当たりの感染者数が最も多い東京と大阪、最も少ない宮城の3都府県で、無作為抽出した20歳以上の男女合わせて7950人を対象に、新型コロナウイルスの抗体検査を実施しました。

その結果、抗体を保有している人の割合は、東京で0.1%、大阪で0.17%、宮城で0.03%だったと16日発表しました。

海外企業のロシュ社とアボット社が製造する2種類の試薬を使った検査で、いずれも陽性となった場合に「陽性」=「抗体を保有している」としました。

https://www3.nhk.or.jp/news/html/20200616/k10012471891000.html

抗体検査についての統計的側面を中心にしたエントリーをあげたのが、4月23日なのでかれこれこの結果を二ヶ月待っていたことになる。

hpo.hatenablog.com

陽性認定の条件が厳しかったのかも知れない。偽陽性を避けるには仕方なかったのだろう。統計的には非常に保守的にデザインされたのだと考える。保守性を加味しても、欧米のような10%以上などの抗体陽性がでなかったことは、日本の新型コロナウィルス対策が結果的に極めて適切であったことが証明された。正直、私は顕在PCR等陽性者1.7万人の50倍から100倍程度は無症状の潜在感染経験者がいるものだと信じていた。

hpo.hatenablog.com

ということで、日本におけるR0が2.5であるという仮説の立証はますます難しくなったと考える。まずはデータの傾向を知るために新規感染者数とRtのグラフを見て、傾向をチェックした。最初に、以前のモデルを0.1%感染経験者という仮定で実データとフィットさせようとしたがどうやってもR0は2以上にはならない。Rtもほとんどの区間で2にはならない。しかも、時系列の先へ行くに従って誤差が大きくなる。1月の新型コロナウィルス上陸から半年近く経っても、この程度の抗体検査陽性、感染経験者数ということは、今回のエピデミックは少なくとも日本においては指数関数的な増加にはなっていない。こう断言する理由も最後で述べる。

f:id:hihi01:20200616224225p:plain

https://docs.google.com/spreadsheets/d/1iChrRlyjpCbJAwwePabR9gIZBfPa-yFCNKrgdgQPYn0/edit?usp=sharing

以前、死亡者数とゴンペルツ曲線の近似を行った。かなりフィットしているので、感染者数もゴンペルツ曲線で近似しうると考えたが、PCR陽性者と潜在感染者の関係がわからなかったので、やらなかった。今回始めて統計的妥当性を持つ調査の結果が明らかになったので、総感染者数の近似を行うことにした。

hpo.hatenablog.com

中国株と欧州株との間で係数に大きな開きがあることは確実なので、屈曲点を仮定して報告美ベースで2月半ばから3月末までの感染者数と5月6日までの感染者数によってゴンペルツ曲線の係数の推定を行った。2月後半から3月末までと、4月1日から5月6日までのそれぞれのデータから係数の同定、推定を行い実データとの比較グラフを作った。当然だが5月初旬まではフィットするがその後は最大感染者数、Gmaxが2万人と過大に推測されたためにあまりフィットしていない。それでも、指数関数、SIR的なモデルの予想よりははるかにフィットしている。指数関数との違いは恣意的にRtあるいはRt=R0 x (1- p)のpを操作しなくとも、時系列の先に行っても実際のデータとの乖離が大きくはならない。

f:id:hihi01:20200616221811p:plain

次に、4月1日前後に存在するであろう屈曲点を求めて、実データ当日と前日の累積感染者数、F(t)とF(t+1)、の縦軸対数グラフにおいて、一次方程式の近似で相関係数R^2が最大となる点を探した。報告日ベースで、3月22日と4月14日当たりでR^2が.99となり、屈曲点でないかと考えた。

まず、4月14日で分けた近似グラフ。黄色の破線は、4月14日から5月6日までのデータで係数を推測している。よくフィットしている。屈曲点の前と後のGmaxが4月14日以前のデータからだと51百万人と、以降のデータからで1.7万人。驚くほどの数字だ。3月の上昇がかなり厳しかったことがうかがえる。ただし、現在1.7万人を超えてしまっているので必ずしても未来予測には使えない。

f:id:hihi01:20200616221440p:plain

次に、同じくR^2が最大.99となった3月22日で同様のことを行った。この時の予想Gmaxは同様に前半データで4,300人、後半データで19百万人。

f:id:hihi01:20200616222405p:plain

この結果を受けて、屈曲点が2つ、同様の傾向を持つ区間が3つ存在すると考察を行った。

この屈曲点をどう考えるか?

欧州株感染者が大量に入国した3月半ばだけが外乱要因であったと考えるとよくつじつまがあう。

押谷先生のインタビュー記事はもう公開が終わってしまったので、画像を再度貼っておく。3月の入国者の感染者数のところは赤線を書き込みたいくらい。

f:id:hihi01:20200617073900p:plain

ちなみに、ゴンペルツ曲線はかなり感染が進んだ後でないと信頼できる係数の算定が難しいという面があるらしい。

よく寝てからまた明日によくよく考える。 → 考えた結果。

ゴンペルツ曲線はある程度の期間データを取らないと近似可能な係数の推定にブレがでるようだし、減衰していった先の方のデータでの一致は指数関数、単純なSIRモデルよりははるかに「まし」に私には見えてるが、それでもだらだらと続く区間については予測が難しいように見える。

K値とゴンペルツ曲線は同じ(形式で表現できる)

本当にきちんと数学を勉強してこなかったことが悔やまれる。こんな明々白々なことがきちんと言葉にできなかった。答えそのものを教えていただいているのに、36時間も目の前に書いてあることが理解できていなかった。恥ずかしい限りだ。

まずは、教えていただいた地下生活者さんのツイットを引用させていただく。

明確に「二重指数関数」と書いてくださっている。もうこの時点で数学の素養のある方には十分だったはず。更に秋山先生のメモを拝見すると、

f:id:hihi01:20200615113531p:plain

http://www.bi.cs.titech.ac.jp/COVID-19/Mathematical_analysis_of_K_indicator.pdf

そして、改めてゴンペルツ曲線を学ばせていただいたBMJrさんの解説を読むと、「ゴンペルツ曲線は二重指数関数」とちゃんと書いてある。対数をきちんと学んでないので直感的にはこのことがわからなかった。

f:id:hihi01:20200615115219p:plain

ゴンペルツ曲線とは何か?(2) - Qiita

従って、この二重指数関数的増加は、初期から飽和状態までに適用できるのは、K値もゴンペルツ曲線でも同様。2u3v4wさんに教えていただいた。

各国に当てはめた中野先生らの論文のグラフにおいても、初期のところの攪乱的な状況と減少していったあとでの「差」が分かる。

f:id:hihi01:20200615115444p:plain

https://www.medrxiv.org/content/10.1101/2020.04.25.20080200v2.full.pdf

BMJrさんもきちんと書いてくださっている。

Gompertz 曲線は、ウシの胎児の子宮内成長曲線だけでなく、全ての生物の成長曲線に適用できる2と言われています。生物の成長曲線の他には、微量の DNA や RNA を増幅する方法(ポリメラーゼ・チェイン・リアクション法:PCR法)による template の増加も S字状を描き、Gompertz 曲線に従うことが報告されています3。PCR法では、サーマルサイクラーの周期に従って指数関数的に template が増加する (exponential growth) はずですが、検査試薬や容器の容量に上限があるため Gompertz 効果(経時的な doubling time の延長からやがてプラトーに達する)が現れるのでしょう。このように「ある倍加時間に伴って指数関数的に増加するが、ある有限の最大値があるために対数表示した時の傾きが時間 t と共にe^−kt倍に減衰する」ものは Gompertz 曲線に従い、累積バグ数にもそのような特性があるのでしょう。

で、次の疑問はK値を考えられた中野先生は生物学におけるゴンペルツ曲線についてご存じなのだろうか?ということになる。しつこいようだが、レヴィット教授は自覚的に当てはめていらっしゃる。中野先生、秋山先生とレヴィット教授の意見交換などを夢視ては失礼だろうか?


■追記

秋山先生がペーパーをアップデートされていた。

なかなか理解が追いついていないがこの辺が肝なのかなと。前後の地下生活者さんの解説をお読みいただきたい。

新型コロナウィルス対策としての健康促進

私は強い経済再生論者だ。先日、健康診断を某所で受けてきた。偶然、安倍首相も人間ドックを受診されたと聞いた。

www.sankei.com

片方で、いまだにあちこちでタバコを吹かしている方々を見かける。若年層でも重篤化の大きなリスク、オッズが指摘されているのに!寛容な私でもびっくりしてしまう。

発表は3月だったが、数少ない重篤化リスクに関する中国の論文でも喫煙は大きなオッズ比を示していた。

f:id:hihi01:20200615082228p:plain

Analysis of factors associated with disease outcomes in hosp... : Chinese Medical Journal

※ 翻訳、記述に基づく表作成は本ブログ。

押谷先生の論文を読んでも結局クラスターの殆どは病院内感染、高齢者施設内感染となっている。

f:id:hihi01:20200615084630p:plain

Early Release - Clusters of Coronavirus Disease in Communities, Japan, January–April 2020 - Volume 26, Number 9—September 2020 - Emerging Infectious Diseases journal - CDC

※ 記述に基づく表作成は本ブログ。

厚労省の大規模調査が待たれるところではあるが、新型コロナウィルスは感染力はそこそこ強いが重篤化、致死率は低い病気なのだ。米国のIHMEによる予測サイトで平気で顕在化しているPCR等陽性者の5倍程度がグラフ上で推測感染者と表記されているのを見ると萎える。

f:id:hihi01:20200615084934p:plain

https://covid19.healthdata.org/japan

日本の現在の状況から言えば、そもそも経済的に第二の緊急事態宣言はほとんど不可能。今度やったら米国、欧州並みの暴動が起こるだろう*1

ということは、「次の波」が来ても、高齢者、因子のある人を隔離して、院内感染、施設内感染を起こさないことが肝要となる。若い方で新型コロナウィルスにかかって亡くなった方は糖尿病だったり、ヘビースモーカーであったと聞く。一般の方々コロナ怖さで病院離れが進んでいるというのは本当に本末転倒。もう国が人間ドック等への補助を徹底して出して、若年層で新型コロナウィルス重篤化因子を持つ人を特定することがどうだろうか?因子保有者は当然、アプリ等で接触者のトラッキングができるようにする。だらには、「次の波」の時はホテル等で自粛することを誓約させる。緊急事態宣言で失われた純粋損失の20兆円、40兆円と比べればその程度の補助など安いものだ。健康診断、人間ドックの全額補助も決してお高くはない。さらには、困窮する病院への支援ともなる。

さらには、タバコの禁煙を呼びかけたいところだが、これは賛否が分かれるだろう。

*1:日本人がというよりも、米国の暴動、事実上の自治権の主張といことは内乱は、第三国の介入を感じる。各国で暴動をしかえる組織、「ダークウェブ」が存在するのではないだろうか?

ゴンペルツ曲線とCovid-19関連死予測

ただいま勉強中。マニアック素人として、ロジスティクス、ロトカ=ヴォルテラ、SIRなどを「いじって」きた者としては定常状態に至る道筋が見えるゴンペルツ曲線は興味深いと。

f:id:hihi01:20200608061750p:plain

ゴンペルツ曲線とは何か?(1) - Qiita

まずは、お手本をそのままスプレッドシートに再現した。

f:id:hihi01:20200608061339p:plain
f:id:hihi01:20200608061547p:plain

ゴンペルツ曲線とは何か?(3) - Qiita

f:id:hihi01:20200608061429p:plain
f:id:hihi01:20200608061619p:plain

Learning Gompertz Curve -> Methematical Prediction of Covid19 death 20200607 - Google スプレッドシート

Gmaxはゴンペルツ曲線の定常状態になった時の予測数。これを日本と米国のCovid19の死者数に当てはめて見た。感染者だと検査数、無症状感染者等の問題があるので。

武漢株(strain)と欧州株の第1波、第2波を分離できなかったので、2月移行の死亡者のデータをそのまま当てはめると、まだ定常状態のはるか手前で最終累計四千人を超える死者がでるまでとなってしまう。感染日ベースで3月末に新規感染者ピークを迎えているので、発症、重篤化、死亡のプロセスは平均で2週間程度だと推測した。そこで、4月11日でデータを分けた。結果、4月11日以降のデータからはほぼ5月後半で定常状態という推測結果に見える。そもそも、新型コロナウィルスによる死亡者がゴンペルツ曲線に合うかどうかからやはり議論しないと確定的なことは言えない。

一方、米国のデータをそのまま推測すると現在10万人超えの死亡者に対して、17万人で定常状態になるという予測となる。感覚以上ではないが、感染曲線、死亡者のデータを眺めていてそのような事態が有り得るようにも見える。こちらも、ゴンペルツ曲線が当てはまるという保証はない。

素人の恐いもの知らずで、この状態でレヴィット教授に質問してしまった。

で、質問してしまったからレヴィット教授のグラフに既に「ゴンペルツ」と書いてあるのに気づいた・・・。 orz....

f:id:hihi01:20200608144406p:plain

数学も、英語も、人口動態ももっと勉強しないとと。


■追記

韓国のMERSの時のデータに今回のスプレッドシートを当てはめるとほぼ、最大感染者数が導けた。

「地下生活者」さんのご指導のおかげ。

https://twitter.com/petenshicomtwitter.com