▶ 学校法人慶應義塾の特許一覧 ▶ 株式会社オトクエストの特許一覧
特開2023-169125抗体価予測モデルおよび抗体価予測モデルを用いたプログラム
(19)【発行国】日本国特許庁(JP)
(12)【公報種別】公開特許公報(A)
(11)【公開番号】P2023169125
(43)【公開日】2023-11-29
(54)【発明の名称】抗体価予測モデルおよび抗体価予測モデルを用いたプログラム
(51)【国際特許分類】
G16H 20/10 20180101AFI20231121BHJP
【ＦＩ】
G16H20/10
【審査請求】未請求
【請求項の数】8
【出願形態】ＯＬ
(21)【出願番号】P 2023080136
(22)【出願日】2023-05-15
(31)【優先権主張番号】P 2022080215
(32)【優先日】2022-05-16
(33)【優先権主張国・地域又は機関】JP
【国等の委託研究の成果に係る記載事項】（出願人による申告）令和３年度、国立研究開発法人日本医療研究開発機構、新興・再興感染症に対する革新的医薬品等開発推進研究事業「多様な前向きコホートを用いたＣＯＶＩＤ－１９ワクチンの多角的解析」委託事業、産業技術力強化法第１７条の適用を受ける特許出願
(71)【出願人】
【識別番号】598121341
【氏名又は名称】慶應義塾
(71)【出願人】
【識別番号】522145177
【氏名又は名称】株式会社オトクエスト
(74)【代理人】
【識別番号】110000062
【氏名又は名称】弁理士法人第一国際特許事務所
(72)【発明者】
【氏名】佐藤 泰憲
(72)【発明者】
【氏名】上蓑 義典
(72)【発明者】
【氏名】長島 健悟
(72)【発明者】
【氏名】村田 満
(72)【発明者】
【氏名】長谷川 直樹
【テーマコード（参考）】
5L099
【Ｆターム（参考）】
5L099AA15
(57)【要約】 （修正有）
【課題】抗体価の経時的変化の予測モデルを構築し、追加接種やそれ以降の接種が必要になる予測時期を提示する抗体価予測モデルおよびそれを用いたプログラムを提供する。
【解決手段】ワクチン接種後の複数の被験者ｉで構成される第１グループに対して所定期間にわたって測定した抗体価の第１検査データに基づいて、ワクチン接種後の抗体価を予測する抗体価予測モデルであって、該抗体価予測モデルは、第１検査データに含まれるデータに関し、被験者ｉのｊ番目の時点をＴｉｊとし、時点Ｔｉｊにおける抗体価をＹｉｊとすると、コンパートメントの容量をａｉｊ、消失速度をｂｉｊとする所定の式に示されるＭコンパートメントモデルを採用し、容量ａｉｊおよび消失速度ｂｉｊに弱情報事前分布を設定し、抗体価Ｙｉｊの対数が他の所定の式に示される対数正規分布に従い、任意のｍコンパートメントの容量ａｉｍおよび消失速度ｂｉｍが多変量正規分布に従う。
【選択図】図１３
【解決手段】ワクチン接種後の複数の被験者ｉで構成される第１グループに対して所定期間にわたって測定した抗体価の第１検査データに基づいて、ワクチン接種後の抗体価を予測する抗体価予測モデルであって、該抗体価予測モデルは、第１検査データに含まれるデータに関し、被験者ｉのｊ番目の時点をＴｉｊとし、時点Ｔｉｊにおける抗体価をＹｉｊとすると、コンパートメントの容量をａｉｊ、消失速度をｂｉｊとする所定の式に示されるＭコンパートメントモデルを採用し、容量ａｉｊおよび消失速度ｂｉｊに弱情報事前分布を設定し、抗体価Ｙｉｊの対数が他の所定の式に示される対数正規分布に従い、任意のｍコンパートメントの容量ａｉｍおよび消失速度ｂｉｍが多変量正規分布に従う。
【選択図】図１３
【特許請求の範囲】
【請求項1】
ワクチン接種後の複数の被験者ｉ（ｉ＝１、２、…、Ｎ）で構成される第１グループに対して所定期間にわたって測定した抗体価の第１検査データに基づいて、ワクチン接種後の抗体価を予測する抗体価予測モデルであって、
前記第１検査データに含まれるデータに関し、被験者ｉのｊ（ｊ＝１、２、…、Ｊ）番目の時点をＴｉｊとし、時点Ｔｉｊにおける抗体価をＹｉｊとすると、
コンパートメントの容量をａｉｊ、消失速度をｂｉｊとする式Ｍ２に示されるＭコンパートメントモデル（Ｍは１以上の整数）を採用し、
前記容量ａｉｊおよび前記消失速度ｂｉｊに弱情報事前分布を設定し、
抗体価Ｙｉｊの対数が式Ｍ１に示される対数正規分布に従い、任意のｍコンパートメントの容量ａｉｍおよび消失速度ｂｉｍが多変量正規分布に従う、ことを特徴とする抗体価予測モデル。
前記第１検査データに含まれるデータに関し、被験者ｉのｊ（ｊ＝１、２、…、Ｊ）番目の時点をＴｉｊとし、時点Ｔｉｊにおける抗体価をＹｉｊとすると、
コンパートメントの容量をａｉｊ、消失速度をｂｉｊとする式Ｍ２に示されるＭコンパートメントモデル（Ｍは１以上の整数）を採用し、
前記容量ａｉｊおよび前記消失速度ｂｉｊに弱情報事前分布を設定し、
抗体価Ｙｉｊの対数が式Ｍ１に示される対数正規分布に従い、任意のｍコンパートメントの容量ａｉｍおよび消失速度ｂｉｍが多変量正規分布に従う、ことを特徴とする抗体価予測モデル。
【数1】
【請求項2】
請求項１に記載の抗体価予測モデルであって、
前記弱情報事前分布は、
前記容量ａｉｊの平均および前記消失速度ｂｉｊの平均は、正規分布に従い、
前記容量ａｉｊの標準偏差および前記消失速度ｂｉｊの標準偏差は、半コーシー分布に従い、
前記容量ａｉｊおよび前記消失速度ｂｉｊの分散共分散行列は、ＬＫＪ分布に従う、ことを特徴とする抗体価予測モデル。
前記弱情報事前分布は、
前記容量ａｉｊの平均および前記消失速度ｂｉｊの平均は、正規分布に従い、
前記容量ａｉｊの標準偏差および前記消失速度ｂｉｊの標準偏差は、半コーシー分布に従い、
前記容量ａｉｊおよび前記消失速度ｂｉｊの分散共分散行列は、ＬＫＪ分布に従う、ことを特徴とする抗体価予測モデル。
【請求項3】
請求項１に記載の抗体価予測モデルであって、
前記被験者ごとの検査データから選別したデータを用いて、非線形回帰モデルを適用して、前記容量ａｉｊおよび前記消失速度ｂｉｊを推定し、
推定された前記容量ａｉｊおよび前記消失速度ｂｉｊのうち、Ｍコンパートメントモデル（式Ｍ２）に当てはめた結果、前記被験者の検査データとほぼ適合する容量ａｉｊおよび前記消失速度ｂｉｊの群を抽出し、
前記抽出された容量ａｉｊおよび消失速度ｂｉｊの群の中でグリッドサーチを行い、初期値となる容量ａおよび消失速度ｂを設定する、ことを特徴とする抗体価予測モデル。
前記被験者ごとの検査データから選別したデータを用いて、非線形回帰モデルを適用して、前記容量ａｉｊおよび前記消失速度ｂｉｊを推定し、
推定された前記容量ａｉｊおよび前記消失速度ｂｉｊのうち、Ｍコンパートメントモデル（式Ｍ２）に当てはめた結果、前記被験者の検査データとほぼ適合する容量ａｉｊおよび前記消失速度ｂｉｊの群を抽出し、
前記抽出された容量ａｉｊおよび消失速度ｂｉｊの群の中でグリッドサーチを行い、初期値となる容量ａおよび消失速度ｂを設定する、ことを特徴とする抗体価予測モデル。
【請求項4】
請求項３に記載の抗体価予測モデルであって、
前記第１検査データを時系列にＫ分割して第１分割データから第Ｋ分割データを形成し、第ｋ（ｋ＝１、…、Ｋ）分割データをテストデータとし、残りの分割データを学習データとし、前記学習データを用いて算出した事後分布と、前記テストデータから任意の時点の抗体価を用い、抗体価の予測される範囲を示す第ｋ予測分布を導出し、前記第ｋ予測分布と前記テストデータとで第ｋ評価をし、
第ｋ評価をｋ＝１からｋ＝Ｋまで行い、
第１評価から第Ｋ評価までに基づく評価をする、ことを特徴とする抗体価予測モデル。
前記第１検査データを時系列にＫ分割して第１分割データから第Ｋ分割データを形成し、第ｋ（ｋ＝１、…、Ｋ）分割データをテストデータとし、残りの分割データを学習データとし、前記学習データを用いて算出した事後分布と、前記テストデータから任意の時点の抗体価を用い、抗体価の予測される範囲を示す第ｋ予測分布を導出し、前記第ｋ予測分布と前記テストデータとで第ｋ評価をし、
第ｋ評価をｋ＝１からｋ＝Ｋまで行い、
第１評価から第Ｋ評価までに基づく評価をする、ことを特徴とする抗体価予測モデル。
【請求項5】
請求項３に記載の抗体価予測モデルであって、
前記第１グループの前記第１検査データに基づいて算出した事後分布を用いて、前記第１グループとは異なる第２グループについて測定した第２検査データの第１時点のデータに基づき、抗体価の予測される範囲を示す予測分布を導出し、
前記予測分布から求めた予測区間に、前記第２検査データに含まれるデータのうち、前記第１時点以降の第２時点のデータが含まれるかどうかを評価する、ことを特徴とする抗体価予測モデル。
前記第１グループの前記第１検査データに基づいて算出した事後分布を用いて、前記第１グループとは異なる第２グループについて測定した第２検査データの第１時点のデータに基づき、抗体価の予測される範囲を示す予測分布を導出し、
前記予測分布から求めた予測区間に、前記第２検査データに含まれるデータのうち、前記第１時点以降の第２時点のデータが含まれるかどうかを評価する、ことを特徴とする抗体価予測モデル。
【請求項6】
請求項１に記載の抗体価予測モデルであって、
前記抗体価予測モデルから、抗体価の予測される範囲を示す予測分布を導出し、
前記予測分布に含まれる予測区間を算出し、
前記予測区間が抗体価の閾値を下回る点を抽出する、ことを特徴とする抗体価予測モデル。
前記抗体価予測モデルから、抗体価の予測される範囲を示す予測分布を導出し、
前記予測分布に含まれる予測区間を算出し、
前記予測区間が抗体価の閾値を下回る点を抽出する、ことを特徴とする抗体価予測モデル。
【請求項7】
請求項６記載の抗体価予測モデルを用いたプログラムであって
ワクチン接種日、抗体価を検査した日、抗体価、抗体検査の方法を入力するステップと、
前記入力したデータから前記抗体価予測モデルによる演算を行うステップと、
前記予測区間が抗体価の前記閾値を下回る点を追加接種やそれ以降の接種の予測時期として出力するステップ、をコンピュータに実行させる、ことを特徴とするプログラム。
ワクチン接種日、抗体価を検査した日、抗体価、抗体検査の方法を入力するステップと、
前記入力したデータから前記抗体価予測モデルによる演算を行うステップと、
前記予測区間が抗体価の前記閾値を下回る点を追加接種やそれ以降の接種の予測時期として出力するステップ、をコンピュータに実行させる、ことを特徴とするプログラム。
【請求項8】
請求項１に記載の抗体価予測モデルであって、
前記ワクチン接種が、コロナウイルスを対象とする、ことを特徴とする抗体価予測モデル。
前記ワクチン接種が、コロナウイルスを対象とする、ことを特徴とする抗体価予測モデル。
【発明の詳細な説明】
【技術分野】
【0001】
本発明は、抗体価予測モデルおよび抗体価予測モデルを用いたプログラムに関する。
【背景技術】
【0002】
ワクチンの接種は、体内にＩｇＭやＩｇＧなどの抗体を産出し、ウイルスに対する感染予防効果および重症化予防効果等をもたらす。従来インフルエンザなどのウイルスに対するワクチン接種は年１回程度の接種にとどまり、接種後の抗体価の推移などを追跡する一般的ニーズもさほどなかった。ところが新型コロナウイルスのパンデミックが長期化するとともに、ワクチン接種後に抗体価が体内で経時的に減少する傾向特性も判明してきたことから、ワクチン接種後に抗体価の推移を予測し、適切な時期に定期的にワクチンの追加接種やそれ以降の接種を行い抗体を増加させる一般的なニーズが公共レベルと個人レベルにおいて高まってきた。
【0003】
また、ワクチンの供給の観点から考えると、ワクチンの生産には製法や製造設備設計が複雑であること、製造管理の難易度が高い等の問題から一般の医薬品に比べて多くの労力と時間がかかる。また、保管方法や有効期間などの点から流通にも慎重を期す必要がある。ワクチン接種をする医師や看護師等の人材の確保も必要である。ワクチン接種の要否を個人が判断するための判断材料としての情報を提供するためには、追加接種やそれ以降の接種が必要になる時期がわかることが望ましい。
【0004】
非特許文献１には、ＣＯＶＩＤ－１９のｍＲＮＡワクチンであるＢＮＴ１６２ｂ２を接種した後２日以内、１４日、２８日、４２日、５６日、および９０日後に採血を行い、１コンパートメント動態モデルを用いて、抗体価が最大濃度になるまでの日数、および抗体価の半減期を推定することが開示されている。
【先行技術文献】
【非特許文献】
【0005】
【非特許文献1】Julien Favresse他、「Antibody titres decline 3-month post-vaccination with BNT 162b2」Emerging Microbes & Infections Volume 10, 2021 - Issue 1
【発明の概要】
【発明が解決しようとする課題】
【0006】
ワクチン接種後の抗体価の経時推移は、年齢や性別などの個人差も含め様々な要因が関係し、本来複雑な推移をたどるものと考えられる。また、それだけに予測モデルをたてるにはできるだけ多くの個人による大量のデータを長期間継続的に取得して、バリデーションを行って精度を高める必要がある。
これに対し非特許文献１で使用している抗体価のモデルは、時間とともに安定的に減衰する単純なモデルである。また３ヶ月以上の長期にわたる抗体価データが十分考慮されておらず、予測モデルの評価（バリデーション）の手法も明らかであるとは言えない。
これに対し非特許文献１で使用している抗体価のモデルは、時間とともに安定的に減衰する単純なモデルである。また３ヶ月以上の長期にわたる抗体価データが十分考慮されておらず、予測モデルの評価（バリデーション）の手法も明らかであるとは言えない。
【0007】
そこで、本発明では、抗体価の経時的変化の予測モデルを構築し、追加接種やそれ以降の接種が必要になる予測時期を提示する技術を提供することを目的とする。
【課題を解決するための手段】
【0008】
上記の課題を解決するために、代表的な本発明の抗体価予測モデルの一つは、ワクチン接種後の複数の被験者ｉ（ｉ＝１、２、…、Ｎ）で構成される第１グループに対して所定期間にわたって測定した抗体価の第１検査データに基づいて、ワクチン接種後の抗体価を予測する抗体価予測モデルであって、前記第１検査データに含まれるデータに関し、被験者ｉのｊ（ｊ＝１、２、…、Ｊ）番目の時点をＴｉｊとし、時点Ｔｉｊにおける抗体価をＹｉｊとすると、コンパートメントの容量をａｉｊ、消失速度をｂｉｊとする式Ｍ２に示されるＭコンパートメントモデル（Ｍは１以上の整数）を採用し、前記容量ａｉｊおよび前記消失速度ｂｉｊに弱情報事前分布を設定し、抗体価Ｙｉｊの対数が式Ｍ１に示される対数正規分布に従い、任意のｍコンパートメントの容量ａｉｍおよび消失速度ｂｉｍが多変量正規分布に従う、ことを特徴とする。
また、追加接種以降の（３回目以降の）接種に対する抗体価予測モデルについては、Ｍコンパートメントモデル（Ｍは１以上の整数）を採用し、それ以外の特徴については上述のモデルと同じとした。
また、追加接種以降の（３回目以降の）接種に対する抗体価予測モデルについては、Ｍコンパートメントモデル（Ｍは１以上の整数）を採用し、それ以外の特徴については上述のモデルと同じとした。
【数1】
【発明の効果】
【0009】
本発明によれば、抗体価の経時的変化の予測モデルを構築し、追加接種やそれ以降の接種が必要になる予測時期を提示することができる。
上記した以外の課題、構成および効果は、以下の実施をするための形態における説明により明らかにされる。
上記した以外の課題、構成および効果は、以下の実施をするための形態における説明により明らかにされる。
【図面の簡単な説明】
【0010】
【発明を実施するための形態】
【0011】
以下、図面を参照して、本発明の実施形態について説明する。なお、この実施形態により本発明が限定されるものではない。また、図面の記載において、同一部分には同一の符号を付して示している。
【0012】
＜第１実施形態＞
（抗体価の検査データ）
図１から図５を参照して、抗体検査の検査結果について説明する。発明者は、所属する機関の職員を検査対象のグループ（以下、「第１グループ」ともいう。）とし、このグループに対してワクチン接種前後に抗体価を測定した。
（抗体価の検査データ）
図１から図５を参照して、抗体検査の検査結果について説明する。発明者は、所属する機関の職員を検査対象のグループ（以下、「第１グループ」ともいう。）とし、このグループに対してワクチン接種前後に抗体価を測定した。
【0013】
図１は、抗体検査を行った週数および実施期間を示す図である。被験者は２０２１年３月にワクチン接種を行い、ワクチン接種から３週後、８週後、１３週後、２６週後に検査を行った。被験者数は６７３人であり、そのうち９５％以上で継続的な検査が行われた。
追加接種後にも検査が実施され、ワクチン追加接種から３週後、８週後、２５週後に検査を行った。
追加接種後にも検査が実施され、ワクチン追加接種から３週後、８週後、２５週後に検査を行った。
【0014】
ワクチンは、Pfizer社／BIONTECH社製の新型コロナウイルスを対象とした感染症の予防のワクチン（SARS-CoV-2）である。抗体の測定に使用した試薬は、Ｓ蛋白ＲＢＤに対するＩｇＧ定量試薬であるAlinity SARS-CoV-2 IgG II Quant、Ｓ蛋白とＮ蛋白それぞれに対するＩｇＧ、ＩｇＭ定量評価のための試薬であるHISCL SARS CoV-2、Pseudo-virusを用いた中和抗体定量評価の試薬（出願人が共同開発）である。
【0015】
抗体価の検査は所定の期間にわたって行われた。まず、ワクチンの接種が行われる前、２０２１年２月２５日から２０２１年３月１２日の間に行われた。ワクチンの接種が２０２１年３月に行われた後、３週目の検査が２０２１年４月１６日から２０２１年４月２８日の間に行われた。続いて、また、ワクチン接種後８週目の検査は、２０２１年５月２０日から２０２１年６月２日までの間に行われた。ワクチン接種後１３週目の検査は、２０２１年５月２０日から２０２１年６月２日までの間に行われた。２６週目の検査は、２０２１年９月２７日から２０２１年１０月８までの間に行われた。
ワクチンの追加接種は２０２１年１２月～２０２２年１月に行われた。ワクチン追加接種後３週目の検査は、２０２２年１月１１日から２０２２年１月２５日の間に行われた。ワクチン追加接種後８週目の検査は、２０２２年２月１０日から２０２２年３月１日までの間に行われた。２５週目の検査は、２０２２年６月２７日から２０２２年７月８までの間に行われた。
ワクチンの追加接種は２０２１年１２月～２０２２年１月に行われた。ワクチン追加接種後３週目の検査は、２０２２年１月１１日から２０２２年１月２５日の間に行われた。ワクチン追加接種後８週目の検査は、２０２２年２月１０日から２０２２年３月１日までの間に行われた。２５週目の検査は、２０２２年６月２７日から２０２２年７月８までの間に行われた。
【0016】
発明者は、第１グループの検査データ（以下、「第１検査データ」という。）を異なる時点ごとに分析した。
【0017】
図２は、検査データを週数ごとの箱ひげ図にした図である。横軸は接種後の週数、縦軸は抗体価の自然対数である。抗体価の単位はＡＵ／ｍＬである。
図３は、ワクチン追加接種後の検査データを週数ごとの箱ひげ図にした図である。横軸は接種後の週数、縦軸は抗体価の自然対数である。抗体価の単位はＡＵ／ｍＬである。
図３は、ワクチン追加接種後の検査データを週数ごとの箱ひげ図にした図である。横軸は接種後の週数、縦軸は抗体価の自然対数である。抗体価の単位はＡＵ／ｍＬである。
【0018】
第３週では、メディアンは１３４６６、第一四分位（Ｑ１）は８４５４、第三四分位（Ｑ３）は２０５１１である。８週目ではメディアンは５４４１、Ｑ１は３６１７、Ｑ３は８３２８である。１３週目ではメディアンは２７４７、Ｑ１は１７８３、Ｑ３は４２５０である。２６週目ではメディアンは８９９、Ｑ１は５７４、Ｑ３は１４６２である。
追加接種後第３週では、メディアンは２０４１９、第一四分位（Ｑ１）は１２５９６、第三四分位（Ｑ３）は３０６５６である。追加接種後８週目ではメディアンは１４５５３、Ｑ１は８９２２、Ｑ３は２２５５７である。追加接種後２５週目ではメディアンは５３８０、Ｑ１は３０５４、Ｑ３は１２２４２である。
追加接種後第３週では、メディアンは２０４１９、第一四分位（Ｑ１）は１２５９６、第三四分位（Ｑ３）は３０６５６である。追加接種後８週目ではメディアンは１４５５３、Ｑ１は８９２２、Ｑ３は２２５５７である。追加接種後２５週目ではメディアンは５３８０、Ｑ１は３０５４、Ｑ３は１２２４２である。
【0019】
どの週数においても、箱で示された四分位範囲の中央近辺に、メディアンが位置している。抗体価のばらつきは非常に少ないことがわかる。なお、２６週目において抗体価が上昇しているデータ群Ｄ１は、被験者にブレイクスルー感染が起こった可能性を示唆していると考えられる。
【0020】
図４は、検査データを年齢と性別ごとに分類したうえで平均をとったものを示した図である。ここでは、対数平均を用いて表示している。年齢は、年齢２３歳から３０歳まで、３０歳から３９歳まで、４０歳から４９歳まで、５０歳から５９歳まで、６０歳から６９歳までの５区分とした。図４Ｍは男性についてのデータであり、図４Ｆは女性についてのデータである。また、閾値５０（対数軸で３．９）は抗体の有無を示す指標であり、例えば抗体価の値が閾値５０を下回る場合には抗体が消失したことを意味する。
【0021】
年齢や性別により値に差が認められるものの、抗体価が減少する傾向は共通していると考えられる。また、大きな傾向をみると、３週目から１３週目までの間の傾きと、１３週目から２６週目の間の傾きでは、減少する傾向が変化している様子が示されている。
【0022】
図５は、検査データを、ワクチン接種後に到達した抗体価の初期値ごとに分類したうえで平均をとったものを示した図である。初期値のデータはワクチン接種後３週目の検査データである。ここでは、検体値の初期値が２００００以上、２００００未満かつ１３０００以上、１３０００未満から８０００以上、８０００未満、の４つに分けたうえで、対数平均を取ったうえで表示している。
【0023】
ワクチンを接種すると、抗体価は一旦、上昇し、その後、減少に転じることがわかる。また、若い人ほど抗体価が高く、女性のほうが男性よりも抗体価の値が高い傾向がある。
抗体価の初期値がいずれの値であっても、抗体価が減少する傾向は共通している。
抗体価の初期値がいずれの値であっても、抗体価が減少する傾向は共通している。
【0024】
上述のような観点から検討をしたのち、発明者は抗体価の推移に関する次のような仮説を立てた。
仮説（１）産生された抗体は指数関数的に減衰し消失する。
なお、ＩｇＧの物理的半減期は構造的な違いによって分類される５つのサブクラス（ＩｇＧ１、ＩｇＧ２ａ、ＩｇＧ２ｂ、ＩｇＧ３、ＩｇＧ４）ごとに異なるが、通常は１５日から２６日程度である。
仮説（２）ｍＲＮＡワクチンの効果で微量の抗体産生が続く可能性もある。常に抗体が足されるため、抗体価の減衰率が低下する。
このことは、例えば、図４および図５において、１３週目の前後で減衰の傾きの変化が、他の週数の前後に見られる変化よりも大きくなっていることから推察される。
仮説（３）追加接種後では、ｍＲＮＡワクチンの効果による微量の抗体産生の反応が定常状態に近づいていると考えられることから減衰率の低下が小さい可能性があり、当てはまりの良いモデルが異なる可能性がある。
仮説（１）産生された抗体は指数関数的に減衰し消失する。
なお、ＩｇＧの物理的半減期は構造的な違いによって分類される５つのサブクラス（ＩｇＧ１、ＩｇＧ２ａ、ＩｇＧ２ｂ、ＩｇＧ３、ＩｇＧ４）ごとに異なるが、通常は１５日から２６日程度である。
仮説（２）ｍＲＮＡワクチンの効果で微量の抗体産生が続く可能性もある。常に抗体が足されるため、抗体価の減衰率が低下する。
このことは、例えば、図４および図５において、１３週目の前後で減衰の傾きの変化が、他の週数の前後に見られる変化よりも大きくなっていることから推察される。
仮説（３）追加接種後では、ｍＲＮＡワクチンの効果による微量の抗体産生の反応が定常状態に近づいていると考えられることから減衰率の低下が小さい可能性があり、当てはまりの良いモデルが異なる可能性がある。
【0025】
（数理モデルの例）
発明者は、抗体価を算出する時点をＴ、抗体価をＹとしたとき、抗体価Ｙの推移を予測するため、予測モデルを検証した。
発明者は、抗体価を算出する時点をＴ、抗体価をＹとしたとき、抗体価Ｙの推移を予測するため、予測モデルを検証した。
【0026】
（コンパートメントモデル）
発明者は、抗体価の減衰メカニズムを解析する過程において、コンパートメントモデルに注目した。
コンパートメントモデルは、薬物が体内に投与されてから、吸収、分布、代謝、排出の体内過程を解析する（ｐｈａｒｍａｃｏｋｉｎｅｔｉｃｓ）ときに用いられるモデルであり、コンパートメントモデルを用いて薬効薬理を予測することが可能となる。
発明者は、抗体価の減衰メカニズムを解析する過程において、コンパートメントモデルに注目した。
コンパートメントモデルは、薬物が体内に投与されてから、吸収、分布、代謝、排出の体内過程を解析する（ｐｈａｒｍａｃｏｋｉｎｅｔｉｃｓ）ときに用いられるモデルであり、コンパートメントモデルを用いて薬効薬理を予測することが可能となる。
【0027】
図６は、薬物が投与される場合の体内過程を示す図である。図６においては体内過程が４つの区間に分けて示されており、体内に投与された薬物は、胃・腸の吸収区間、血液による分布区間、肝臓による代謝区間、腎臓による排泄区間を経て、体外に排泄される。このような体内過程を解析するにあたり、コンパートメントモデルは、「ある区間から薬物が排出される速度は、その区間における薬物濃度に比例する」という物理的仮定を設定する。これを一次速度過程という。なお、薬物の移行速度が時間によって変化せず一定である場合は、０次速度過程という。
【0028】
図７は、静注１コンパートメントモデルを示す図である。静注１コンパートメントモデルは、薬物が静脈内に投与された場合の体内過程を解析するために用いられるモデルであり、例えば、血液中の薬と組織中の薬が速やかに濃度平衡に到達するような場合に、体内の薬物量を近似するのに用いられる。ここで、薬物の血中濃度をＣ（ｔ）とすると、薬物の排出速度はｄＣ（ｔ）／ｄｔと表すことができる。排出速度に関する定数をｋとおくと、静注１コンパートメントモデルの一次速度過程について次の式１が導かれる。
【0029】
【数2】
【0030】
【0031】
図８に示されるような指数関数的な減衰の様子が示される場合には、静注１コンパートメントモデルを適用することができると考えられる。例えば、検査データとして、薬物の投与された量をｑ０、血中濃度Ｃがある場合には、静注１コンパートメントモデルを適用して、排出速度Ｖを推定することができる。
【0032】
（実施例１…２（Ｍ）コンパートメント型の階層ベイズモデル）
上述の仮説（１）は、抗体価が指数関数的に減衰するとした。仮説（２）では、減衰率が経時的に変化するとした。減衰率の変化には、複数のメカニズムが影響する可能性も考えられる。静注１コンパートメントモデルでは、これらの減衰率の変化を十分に考慮することができない問題がある。
上述の仮説（１）は、抗体価が指数関数的に減衰するとした。仮説（２）では、減衰率が経時的に変化するとした。減衰率の変化には、複数のメカニズムが影響する可能性も考えられる。静注１コンパートメントモデルでは、これらの減衰率の変化を十分に考慮することができない問題がある。
【0033】
そこで、発明者は、Ｍを１以上の整数とし、抗体価の減少の様子にＭコンパートメントモデルを仮定した。減衰の様子が変化する点や、抗体産生が続く可能性がある点を、２以上のコンパートメントを用いて表現することを試みた。
【0034】
ここで、Ｍ＝２とした２コンパートメントモデルを用いる場合を説明する。２コンパートメント型の階層ベイズモデルは、抗体価の経時変化を２コンパートメントモデル（式２ｃ－１）に近似させた場合、抗体価Ｙが平均を２コンパートメントモデルｆ２（Ｔ）、標準偏差σＹとした対数正規分布（式２ｃ－２）に従うとしたものである。２コンパートメントの場合の初期分布と減衰速度のパラメータ（ａ、ｂ）（ｃ、ｄ）は、多変量正規分布に従うと考える（式２ｃ－３）。コンパートメントの容量ａおよびｃは抗体価の初期値に関連するパラメータであり、排出速度に相当するｂおよびｄは抗体価の消失速度に関連するパラメータである。ａとｂには相関があると考えて、ａとｂそれぞれの第１グループ内の平均（μａとμｂ）は多変量正規分布に従うとした（式２ｃ－３）。ｃとｄについてもａとｂと同様に多変量正規分布に従うとした（式２ｃ－４）。事後分布を求める際に、事後分布に大きな影響を与えない弱情報事前分布を与える事が望ましい。個々人のデータに非線形回帰モデルを当てはめた結果から予想される値よりも十分に大きな値となるような弱情報事前分布を考慮して、ａとｂの全体平均は平均０で標準偏差１００の正規分布に従うとし、ａとｂの標準偏差が位置パラメータ０、尺度パラメータ５０の半コーシー分布に従うとする。分散共分散行列Σには、ＬＫＪｃｏｒｒ＿Ｃｈｏｌｅｓｋｙ分布（以下、「ＬＫＪ分布」ともいう。）に従うとし、形状パラメータを１として、無情報の一様分布に相当するものとした（式２ｃ－５）。ｃとｄについても、ａとｂと同様の弱情報事前分布を設定した（式２ｃ－６）。なお、式２ｃ－２などの「Ｎｏｒｍａｌ（）」および式２ｃ－５の「Ｎ（）」は正規分布、式２ｃ－３など「Ｍｕｌｔｉｎｏｒｍａｌ（）」は多変量正規分布、式２ｃ－５などの「ｈａｌｆｃａｕｃｈｙ（）」は半コーシー分布、をそれぞれ示す。
【0035】
【数3】
【0036】
（実施例２…ワイブル型の階層ベイズモデル）
仮説（１）および仮説（２）を表現できる可能性がある分布として、発明者はワイブル分布に着目した。ワイブル分布は製品寿命を示す場合に用いられることが多く、指数関数の項を含んでいる。ワイブル型の階層ベイズモデルは、抗体価Ｙの経時変化をワイブル分布（式ｗ－１）に近似させた場合、抗体価Ｙが平均をワイブル分布とする対数正規分布に従うとする。パラメータｂとｄに関して、ｄは抗体価の初期値に相当し、ｂは消失速度に相当するものであり、１コンパートメント型階層ベイズモデルの場合と同じである。
仮説（１）および仮説（２）を表現できる可能性がある分布として、発明者はワイブル分布に着目した。ワイブル分布は製品寿命を示す場合に用いられることが多く、指数関数の項を含んでいる。ワイブル型の階層ベイズモデルは、抗体価Ｙの経時変化をワイブル分布（式ｗ－１）に近似させた場合、抗体価Ｙが平均をワイブル分布とする対数正規分布に従うとする。パラメータｂとｄに関して、ｄは抗体価の初期値に相当し、ｂは消失速度に相当するものであり、１コンパートメント型階層ベイズモデルの場合と同じである。
【0037】
【数4】
【0038】
（実施例３…二重指数関数型の階層ベイズモデル）
仮説（１）および仮説（２）を表現できる可能性がある分布として、発明者は二重指数関数型の階層ベイズモデルに着目した。二重指数関数は、指数関数の肩に指数関数を持つ。二重指数関数型の階層ベイズモデルは、抗体価Ｙの経時変化を二重指数関数（式ｄ－１）に近似させた場合、抗体価Ｙが平均が二重指数関数型とする対数正規分布に従うとする。
仮説（１）および仮説（２）を表現できる可能性がある分布として、発明者は二重指数関数型の階層ベイズモデルに着目した。二重指数関数は、指数関数の肩に指数関数を持つ。二重指数関数型の階層ベイズモデルは、抗体価Ｙの経時変化を二重指数関数（式ｄ－１）に近似させた場合、抗体価Ｙが平均が二重指数関数型とする対数正規分布に従うとする。
【0039】
【数5】
【0040】
（実施例４…１コンパートメント型の階層ベイズモデル）
追加接種後については、仮説（３）から減衰率の低下が小さい可能性がある。仮説（３）を表現できる可能性がある分布として、発明者は１コンパートメント型の階層ベイズモデルに着目した。１コンパートメント型の階層ベイズモデルは、減衰率の低下が無視できるほど小さい場合と考えることができる。１コンパートメント型の階層ベイズモデルは、抗体価Ｙの経時変化が１コンパートメントモデル（式１ｃ－１）に近似させた場合、抗体価Ｙが平均を１コンパートメントモデルｆ１（Ｔ）、標準偏差σＹとした対数正規分布（式１ｃ－２）に従うとする。検査データのグループ差や個人差による変動分が正規分布に従うと考えている。２コンパートメント型の階層ベイズモデルのうち１コンパートメント分に該当するものであり、事前分布（式１ｃ－４）は２コンパートメント型の階層ベイズモデルの場合と同じである。
追加接種後については、仮説（３）から減衰率の低下が小さい可能性がある。仮説（３）を表現できる可能性がある分布として、発明者は１コンパートメント型の階層ベイズモデルに着目した。１コンパートメント型の階層ベイズモデルは、減衰率の低下が無視できるほど小さい場合と考えることができる。１コンパートメント型の階層ベイズモデルは、抗体価Ｙの経時変化が１コンパートメントモデル（式１ｃ－１）に近似させた場合、抗体価Ｙが平均を１コンパートメントモデルｆ１（Ｔ）、標準偏差σＹとした対数正規分布（式１ｃ－２）に従うとする。検査データのグループ差や個人差による変動分が正規分布に従うと考えている。２コンパートメント型の階層ベイズモデルのうち１コンパートメント分に該当するものであり、事前分布（式１ｃ－４）は２コンパートメント型の階層ベイズモデルの場合と同じである。
【0041】
【数6】
【0042】
図９は、各モデルを適用した場合の抗体価Ｙの経時変化を模式的に示す図である。図９は、式２ｃ－１、式ｗ－１、式ｄ－１、式１ｃ－１を模式的に示している。図９（ａ）は経時的に変化する様子を示し、図９（ｂ）は対数で経時変化を示す。凡例２－ｃｏｍｐは２コンパートメント型の階層ベイズモデルを適用した場合であり、凡例Ｗｅｉｂｕｌｌはワイブル型の階層ベイズモデルを適用した場合であり、凡例ｄｌｂ．ｅｘｐは二重指数関数型の階層ベイズモデル適用した場合であり、凡例１－ｃｏｍｐは１コンパートメント型の階層ベイズモデルを適用した場合である。
【0043】
（事後分布の導出方法…Ｍコンパートメントの階層ベイズモデル）
ここで、一般的なＭコンパートメントの階層ベイズモデルを用いて、事後分布を導出する方法を説明する。
ここで、一般的なＭコンパートメントの階層ベイズモデルを用いて、事後分布を導出する方法を説明する。
【0044】
【数7】
【0045】
式Ｍ１および式Ｍ２は、ワクチン接種後の複数の被験者ｉ（ｉ＝１、２、…、Ｎ）で構成される第１グループに対して所定期間にわたって測定した抗体価の第１検査データに基づいて、ワクチン接種後の抗体価を予測する抗体価予測モデルである。ここで、第１検査データに含まれるデータに関し、被験者ｉのｊ（ｊ＝１、…、Ｊ）番目の時点をＴｉｊとし、点Ｔｉｊにおける抗体価をＹｉｊとする。コンパートメントの容量をａｉｊ、消失速度をｂｉｊとする式Ｍ２に示されるＭコンパートメントモデル（Ｍは１以上の整数）を採用し、抗体価Ｙｉｊの対数が式Ｍ１に示される対数正規分布に従い、任意のｍコンパートメントの容量ａｉｍおよび消失速度ｂｉｍが多変量正規分布（式Ｍ７）に従う、ことを仮定する。式Ｍ２における指示変数ｌｍはコンパートメント数を定めるパラメータである。ｍコンパートメントモデル（ｍは１以上の整数とする。）を用いる場合は、Ｉ１＝Ｉ２＝…＝Ｉｍ＝１とし、Ｉｍ＋１＝…＝ＩＭ＝０とする。なお、式Ｍ４などに用いられている記号「´」はベクトルまたは行列の転置、式Ｍ１３の「ｄｉａｇ」は対角行列、式Ｍ１１の「ｖｅｃｈ（）」はｖｅｃｈ作用素、をそれぞれ示している。
【0046】
ｍコンパートメントの容量（ａｉｍ）（式Ｍ６）と消失速度（ｂｉｍ）（式Ｍ７）をまとめたパラメータベクトルθｉ（式Ｍ４）は被験者ごとに異なる値を持ち、θｉの値により各個人の抗体価の平均的な推移が規定される。また、ｍ番目のコンパートメントのパラメータの組（ａｉｍ、ｂｉｍ）は、平均ベクトルがμｍで分散共分散行列がΣｍの多変量正規分布に従う階層構造を持つことを仮定した（式Ｍ７）。
【0047】
さらに、μｍとΣｍに与える事前分布の超パラメータについては、事後分布に大きな影響を与えない弱情報の事前分布を与える事が望ましく、被験者ごとのデータに非線形回帰モデルを当てはめた結果から予想されるθの取る範囲よりも十分に広くなるような超パラメータとして、次の式Ｍ１４から式Ｍ１８を仮定した。ここでは、容量ａｉｊおよび消失速度ｂｉｊに弱情報事前分布を設定したことを示している。式Ｍ１４および式Ｍ１５に示されるように容量ａｉｊの平均μａｍおよび消失速度ｂｉｊの平均はμｂｍは正規分布に従い、式Ｍ１６および式Ｍ１７に示されるように容量ａｉｊの標準偏差σａｍおよび消失速度ｂｉｊの標準偏差σｂｍは半コーシー分布に従い、式Ｍ１８に示されるように容量ａｉｊおよび消失速度ｂｉｊの分散共分散行列ΣｍはＬＫＪ分布に従う。
【0048】
【数8】
【0049】
上述のモデルにおいて、同時分布の対数をとったもの（以下、「尤度関数」ともいう。）Ｌ（θ）を式Ｍ１９に示す。ここで、右辺第１項pＹはＹｉｊの確率密度関数、右辺第２項pａｂは（ａｉｍ、ｂｉｍ）の同時確率密度関数、右辺第３項（pμａ、pμｂ、pσａ、pσｂ、pΣ）それぞれは事前分布の確率密度関数である。また、∇θＬ（θ）はＬ（θ）をベクトルθで微分したものとし、モデルにおけるパラメータθの要素数をνとする。
【0050】
【数9】
【0051】
【0052】
ステップＳ１において、パラメータθに対する初期値θ０と、微小な定数ε、定数Ｓ、Ｕを与える。
ステップＳ２において、状態ｓ＝１を設定する。
ステップＳ３において、平均ベクトルが０で分散共分散行列が単位行列Ｉνのν次元多変量正規分布に従う乱数を生成し、ｒ０とする。
ステップＳ４において、θｓ＝θｓ－１、θ(～)＝θｓ－１、ｒ(～)＝ｒ０を設定する。
ステップＳ２において、状態ｓ＝１を設定する。
ステップＳ３において、平均ベクトルが０で分散共分散行列が単位行列Ｉνのν次元多変量正規分布に従う乱数を生成し、ｒ０とする。
ステップＳ４において、θｓ＝θｓ－１、θ(～)＝θｓ－１、ｒ(～)＝ｒ０を設定する。
【0053】
ステップＳ５において、状態ｕ＝１を設定する。
ステップＳ６において、ｒ(～)´＝ｒ(～)＋（ε／２）∇θＬ（θ(～)）、θ(～)ｕ＝θ(～)＋εｒ(～)´、ｒ(～)ｕ＝ｒ(～)´＋（ε／２）∇θＬ（θ(～)ｕ）を求める。
ステップＳ６において、ｒ(～)´＝ｒ(～)＋（ε／２）∇θＬ（θ(～)）、θ(～)ｕ＝θ(～)＋εｒ(～)´、ｒ(～)ｕ＝ｒ(～)´＋（ε／２）∇θＬ（θ(～)ｕ）を求める。
【0054】
ステップＳ７において、ｕにｕ＋１を代入する。
ステップＳ８において、状態ｕがＵより大きいかどうか判定する。状態ｕがＵ以下である場合、ステップＳ５に戻る。状態ｕがＵより大きい場合、ステップＳ９に移行する。
ステップＳ９において、確率α＝｛１、ｅｘｐ｛Ｌ（θ(～)Ｕ）＋（１／２）（ｒ(～)Ｕ´）ｒ(～)Ｕ｝／ｅｘｐ｛Ｌ（θｓ－１）＋（１／２）（ｒ０）´ｒ０｝｝で、θｓ＝θＵおよびｒｓ＝ｒ(～)Ｕに更新する。
ステップＳ１０において、状態ｓがＳより大きいかどうか判定する。状態ｓがＳ以下である場合、ステップＳ２に戻る。状態ｓがＳより大きい場合、アルゴリズムを終了する。
ステップＳ８において、状態ｕがＵより大きいかどうか判定する。状態ｕがＵ以下である場合、ステップＳ５に戻る。状態ｕがＵより大きい場合、ステップＳ９に移行する。
ステップＳ９において、確率α＝｛１、ｅｘｐ｛Ｌ（θ(～)Ｕ）＋（１／２）（ｒ(～)Ｕ´）ｒ(～)Ｕ｝／ｅｘｐ｛Ｌ（θｓ－１）＋（１／２）（ｒ０）´ｒ０｝｝で、θｓ＝θＵおよびｒｓ＝ｒ(～)Ｕに更新する。
ステップＳ１０において、状態ｓがＳより大きいかどうか判定する。状態ｓがＳ以下である場合、ステップＳ２に戻る。状態ｓがＳより大きい場合、アルゴリズムを終了する。
【0055】
上述のアルゴリズムにおいて、状態ｓはハミルトニアン・モンテカルロ法のステップ数を示している。ステップが進むごとに、θ（式Ｍ１０）の事後分布から、θの１組がサンプリングされる。また、ステップＳ５からステップＳ８は、リープフロッグ積分を用いて微分方程式を数値積分する工程を示しており、状態ｕはリープフロッグ積分におけるステップ数を表す変数である。ステップＳ５からステップＳ８を繰り返すことで、初期状態（ｕ＝０）から状態εＵまで進んだ仮想的な時点でのθとｒを求めることができる。
【0056】
（初期値の設定方法）
ステップＳ１の初期値θ０の設定について、θの初期値はＬ（θ）や∇θＬ（θ）の算出に影響を及ぼすと考えられる。例えば、Ｌ（θ）の密度が集中していない値、Ｌ（θ）がオーバーフローするような値、勾配∇θＬ（θ）が消失するような値をθの初期値として与えると、アルゴリズムが実際の事後分布とずれた場所でθの抽出をするか、計算不能になるか、同じ場所を移動することになり、事後分布に適切に収束しない恐れがある。このような値を避けて初期値を設定する必要がある。
ステップＳ１の初期値θ０の設定について、θの初期値はＬ（θ）や∇θＬ（θ）の算出に影響を及ぼすと考えられる。例えば、Ｌ（θ）の密度が集中していない値、Ｌ（θ）がオーバーフローするような値、勾配∇θＬ（θ）が消失するような値をθの初期値として与えると、アルゴリズムが実際の事後分布とずれた場所でθの抽出をするか、計算不能になるか、同じ場所を移動することになり、事後分布に適切に収束しない恐れがある。このような値を避けて初期値を設定する必要がある。
【0057】
実データから外れた値を入れるとＬ（θ）や∇θＬ（θ）の値は実際の事後分布とずれると考えられるため、まずは個々の被験者にコンパートメントモデルを当てはめてパラメータを確定した結果を得る。実際のデータには欠測などがあることから、必ずしも個々の被験者のパラメータをうまく推定できるわけではないので、個々の被験者のパラメータとしてうまく推定できた値だけを抽出し、θの候補範囲を絞り込む。この絞り込んだ候補範囲に対してグリッドサーチを行い、色々な値を代入して得たＬ（θ）や∇θＬ（θ）の数値を選別して（オーバーフローしそうな値や異常値を除去する等）、最終的なθの初期値の範囲を決定する。
【0058】
図１１は、事後分布を求めるときに用いる初期値の設定方法を示すフローチャートである。
ステップＳ１１において、被験者ごとの検査データから選別したデータを用いて、容量ａｉおよび消失速度ｂｉを推定する。具体的には、被験者ｉの３週目の抗体価、８週目の抗体価、１３週目の抗体価、２６週目の抗体価の散布図を描き、この散布図にＭコンパートメントモデル（式Ｍ２）の非線形回帰モデルを適用して、容量ａｉｍと消失速度ｂｉｍを推定する。また、追加接種後については、被験者ｉの３週目の抗体価、８週目の抗体価、２５週目の抗体価の散布図から、同様の方法により初期値を設定する。
ステップＳ１１において、被験者ごとの検査データから選別したデータを用いて、容量ａｉおよび消失速度ｂｉを推定する。具体的には、被験者ｉの３週目の抗体価、８週目の抗体価、１３週目の抗体価、２６週目の抗体価の散布図を描き、この散布図にＭコンパートメントモデル（式Ｍ２）の非線形回帰モデルを適用して、容量ａｉｍと消失速度ｂｉｍを推定する。また、追加接種後については、被験者ｉの３週目の抗体価、８週目の抗体価、２５週目の抗体価の散布図から、同様の方法により初期値を設定する。
【0059】
ステップＳ１２において、初期値となりうる容量および消失速度の候補を絞り込む。具体的には、ステップＳ１１において推定された容量ａｉおよび消失速度ｂｉのうち、Ｍコンパートメントモデル（式Ｍ２）に当てはめた結果、被験者の検査データとほぼ適合する容量ａｉおよび消失速度ｂｉの群を抽出する。抽出の方法としては、例えば、四分位範囲（ＩＱＲ）を用いることが可能である。算出された容量と消失速度を順に並べて、第１四分位数と第３四分位数を算出し、四分位範囲（ＩＱＲ）＝第３四分位数―第１四分位数とする。外れ値＜第１四分位数―（１．５×ＩＱＲ）または外れ値＞第３四分位数＋（１．５×ＩＱＲ）とし、外れ値となった容量と消失速度を除いたものを抽出する。
【0060】
ステップＳ１３において、抽出された容量ａｉおよび消失速度ｂｉの群の中でグリッドサーチをする。抽出された容量と消失速度の組合せをＬ（θ）や∇θＬ（θ）に代入し、Ｌ（θ）や∇θＬ（θ）がオーバーフローしそうな値や異常値を取らない適切な値となる容量と消失速度をθの初期値の範囲として設定する。設定した初期値の範囲から抽出した疑似乱数の組を初期値θ０として設定する。
【0061】
なお、初期値の設定の方法は図１１および上述の方法に限定されない。事後分布が収束するように、初期値の設定の方法を適宜選択することができる。
【0062】
（４モデルの適用結果）
図１２は、４モデルの情報量基準（ＬＯＯＩＣ： Ｌｅａｖｅ－Ｏｎｅ－Ｏｕｔ Ｉｎｆｏｒｍａｔｉｏｎ Ｃｒｉｔｅｒｉｏｎ）を示す図である。情報量基準の値が最も小さい２コンパートメント型の階層ベイズモデルが、４モデルのうちで適合度がもっともよい。この理由としては、薬物の場合では血中濃度が急激に減少するα相と徐々に血中濃度が下がっていくβ相があるところ、抗体価の減少についても図４および図５に１３週目の前後で２相ごとの消失の違いをとらえたモデリングが実現できたためと考えられる。このことからは、発明者は、抗体価推移には、Ｍが１以上のＭコンパートメント型の階層ベイズモデルを用いることが適切であると考える。また、１コンパートメント型よりも二重指数型の方が適合度がよい理由は、減衰率が低下していく様子をより適切に考慮できたモデルであったためと考えられる。
図１２には、追加接種後の４モデルのＬＯＯＩＣも示している。追加接種後については、情報量基準の値が最も小さい１コンパートメント型の階層ベイズモデルが、４モデルのうちで適合度がもっともよい。この理由としては、追加接種後では、複数回のｍＲＮＡワクチン接種によって、微量の抗体産生の反応が定常状態に近づき、減衰率の低下が非常に小さく無視できる大きさとなっており、一定の消失をとらえたモデリングが実現できたためと考えられる。このことからは、発明者は、追加接種後の抗体価推移には、Ｍが１以上のＭコンパートメント型の階層ベイズモデルを用いることが適切であると考える。
図１２は、４モデルの情報量基準（ＬＯＯＩＣ： Ｌｅａｖｅ－Ｏｎｅ－Ｏｕｔ Ｉｎｆｏｒｍａｔｉｏｎ Ｃｒｉｔｅｒｉｏｎ）を示す図である。情報量基準の値が最も小さい２コンパートメント型の階層ベイズモデルが、４モデルのうちで適合度がもっともよい。この理由としては、薬物の場合では血中濃度が急激に減少するα相と徐々に血中濃度が下がっていくβ相があるところ、抗体価の減少についても図４および図５に１３週目の前後で２相ごとの消失の違いをとらえたモデリングが実現できたためと考えられる。このことからは、発明者は、抗体価推移には、Ｍが１以上のＭコンパートメント型の階層ベイズモデルを用いることが適切であると考える。また、１コンパートメント型よりも二重指数型の方が適合度がよい理由は、減衰率が低下していく様子をより適切に考慮できたモデルであったためと考えられる。
図１２には、追加接種後の４モデルのＬＯＯＩＣも示している。追加接種後については、情報量基準の値が最も小さい１コンパートメント型の階層ベイズモデルが、４モデルのうちで適合度がもっともよい。この理由としては、追加接種後では、複数回のｍＲＮＡワクチン接種によって、微量の抗体産生の反応が定常状態に近づき、減衰率の低下が非常に小さく無視できる大きさとなっており、一定の消失をとらえたモデリングが実現できたためと考えられる。このことからは、発明者は、追加接種後の抗体価推移には、Ｍが１以上のＭコンパートメント型の階層ベイズモデルを用いることが適切であると考える。
【0063】
（予測分布の計算）
図１３は、２コンパートメント型の階層ベイズモデルから導出された予測分布の第１例を示す図である。バンドｐｉ１は、抗体価の予測される範囲を示す予測分布に従う値のうち２．５％点から９７．５に該当する範囲であり、９５％予測区間を示す。グラフｇ１は、抗体価が閾値５０を切る付近を拡大したものである。
図１３は、２コンパートメント型の階層ベイズモデルから導出された予測分布の第１例を示す図である。バンドｐｉ１は、抗体価の予測される範囲を示す予測分布に従う値のうち２．５％点から９７．５に該当する範囲であり、９５％予測区間を示す。グラフｇ１は、抗体価が閾値５０を切る付近を拡大したものである。
【0064】
図１３において、検査データに含まれる３週目、８週目、１３週目、２６週目の値が実測値である。２６週目以降のデータは予測値である。バンドｐｉ１は４６週目のあたりで閾値５０を下回る。図１３の検査データの被験者は、４６週のあたりで抗体価が閾値５０を下回り、抗体が消失することが予測できる。
【0065】
図１４は、２コンパートメント型の階層ベイズモデルから導出された予測分布の第２例を示す図である。図１４は抗体価が初期値で５００００程度の値を取った被験者の検査データを用いたものである。３週目、８週目、１３週目、２６週目の値が実測値であり、２６週目以降のデータは予測値である。バンドｐｉ２は９５％予測区間を示す。図１４の検査データの被験者は、９０週のあたりで抗体価が閾値５０を下回ることが予測できる。
【0066】
図１５は、２コンパートメント型の階層ベイズモデルから導出された予測分布の第３例を示す図である。図１５はデータが欠測しており３週目のみ実測値がある場合の検査データを用いたものである。バンドｐｉ３は９５％予測区間を占めず。図１５の検査データの被験者は、４０週のあたりで抗体価が閾値５０を下回ることが予測できる。
【0067】
なお、検査対象としたグループの特性、被験者の特性のほか、検査方法によっても抗体価に違いが生じる場合がある。この場合には、所定の関係式を用いて抗体価を補正したうえで予測モデルによる演算をすることによって、適切な予測値を得ることができる。
【0068】
（予測モデルのバリデーション）
図１６を用いて、予測モデルのバリデーションについて説明する。図１６は、バリデーションの手法を模式的に示す図である。バリデーションは、内的妥当性の検証と外的妥当性の検証である。
図１６を用いて、予測モデルのバリデーションについて説明する。図１６は、バリデーションの手法を模式的に示す図である。バリデーションは、内的妥当性の検証と外的妥当性の検証である。
【0069】
内的妥当性の検証は、クロスバリデーション（交差検証）によって行う。図１６（ａ）は、１０分割のクロスバリデーションを説明する図である。まず、第１検査データを時系列に１０分割して第１分割データから第１０分割データを形成する。試行１では第１分割データをテストデータとし、残りの第２分割データから第１０分割データまでをテストデータとする。学習用データを用いて算出した事後分布と、テストデータに含まれる被験者のデータから、３週、８週、１３週の抗体価をランダムに１点だけ用い、第１予測分布を導出する。そして、第１予測分布とテストデータとで後述の評価を行う。残りの分割データそれぞれについてもテストデータとした試行２から試行１０を行い、結果として導出した予測分布を評価対象とした。評価は、評価（１）として２６週時点の予測分布における９５％予測区間がテストデータの実測値を含む、評価（２）として予測分布における９５％予測区間の下限値が＜閾値５０となる時点においてテストデータの実測値も＜閾値５０となる、の２観点で行った。結果として、評価（１）については、９７．０％のテストデータが評価観点に該当した。評価（２）については、９９．５％のテストデータが評価観点に該当した。
また、追加接種後についても同様の方法により、内的妥当性の検証を実施した。結果として、評価（１）については、９３．８％のテストデータが評価観点に該当した。評価（２）については、５０を下回った実測値が存在しなかったため、評価できなかった。
また、追加接種後についても同様の方法により、内的妥当性の検証を実施した。結果として、評価（１）については、９３．８％のテストデータが評価観点に該当した。評価（２）については、５０を下回った実測値が存在しなかったため、評価できなかった。
【0070】
外的妥当性の検証は、第１検査データから算出された事後分布を用いて導出された予測分布の範囲内に、外部データが含まれるかどうかで行った。外部データとして、発明者が所属する機関とは異なる機関に所属する被検者グループについて測定した検査データ（以下、「第２検査データ」という。）を用いた。図１６（ｂ）は、外部妥当性の検証を説明する図である。
ここで、第１検査データと第２検査データとでは測定キットが異なっていた。第２検査データに用いたＢＡＵ／ｍＬ単位は、適当な係数を乗算してＡＵ／ｍＬ単位に変換した。第１検査データによって構築した予測モデルに、第２検査データの３週目～４週目の抗体価を入力して、第２検査データに基づく予測分布を導出した。評価は、１２週時点の９５％予測区間が実測値を含む、の観点で行った。結果として、９２．３％の検査データで評価観点に該当した。
ここで、第１検査データと第２検査データとでは測定キットが異なっていた。第２検査データに用いたＢＡＵ／ｍＬ単位は、適当な係数を乗算してＡＵ／ｍＬ単位に変換した。第１検査データによって構築した予測モデルに、第２検査データの３週目～４週目の抗体価を入力して、第２検査データに基づく予測分布を導出した。評価は、１２週時点の９５％予測区間が実測値を含む、の観点で行った。結果として、９２．３％の検査データで評価観点に該当した。
【0071】
（作用・効果）
本発明のＭコンパートメント型の階層ベイズモデルは、被験者数６７３人の検査データに基づいて予測分布を導出している。このような膨大な検査データに基づく予測分布であるため、個人差等の影響も十分に反映されていると考えられる。またモデルが収束するようデータの選別（前処理）も行われた。こうしたことから他のモデルと比較して検査データへの適合度もよく、内部妥当性および外部妥当性も十分に担保されており、予測精度は大きく改善されている。したがって、本発明を抗体価予測モデルに適用することによって、抗体価の経時的変化を精度よく表すことができ、追加接種やそれ以降の接種が必要になる予測時期を得ることができる。
本発明のＭコンパートメント型の階層ベイズモデルは、被験者数６７３人の検査データに基づいて予測分布を導出している。このような膨大な検査データに基づく予測分布であるため、個人差等の影響も十分に反映されていると考えられる。またモデルが収束するようデータの選別（前処理）も行われた。こうしたことから他のモデルと比較して検査データへの適合度もよく、内部妥当性および外部妥当性も十分に担保されており、予測精度は大きく改善されている。したがって、本発明を抗体価予測モデルに適用することによって、抗体価の経時的変化を精度よく表すことができ、追加接種やそれ以降の接種が必要になる予測時期を得ることができる。
【0072】
【0073】
ステップＳ２１において、検査データを収集する。検査データには、抗体価のほか、被験者を識別する情報、年齢、性別、検査日時、検査方法、などが含まれる。
【0074】
ステップＳ２２において、検査データから抗体価の減衰の特性を判定する。例えば、抗体価の対数平均をとったうえで、１回目の検査と２回目の検査の間の抗体値の減少割合ｒ１、２回目の検査と３回目の検査の間の抗体値の減少割合ｒ２・・・、ｎ回目の検査と（ｎ＋１）回目の検査の間の抗体値の減少割合ｒｎ＋１について分析する。上述の分析は、年齢、性別、抗体価の初期値、等の観点に分けて行ってもよい。
【0075】
ステップＳ２３において、Ｍコンパートメント型の階層ベイズモデルを設定する。ここで、Ｍの値は２に限定されない。例えば、ステップＳ２２における分析に基づいて抗体価減少がどのような変化を持つか判定し、判定の結果にもとづいてＭの値を設定してもよい。
【0076】
ステップＳ２４において、Ｍコンパートメントモデルのパラメータの初期値を設定する。設定方法としては、図１１に示す方法を適用することが可能である。
【0077】
ステップＳ２５において、Ｍコンパートメントモデルのパラメータの事後分布を算出する。第１実施形態においては、リープフロッグ積分を用いてハミルトニアン・モンテカルロ法により事後分布を算出したが、算出方法はこれに限定されない。
【0078】
ステップＳ２６において、追加接種またはそれ以降の接種予定者のデータを入力する。入力されるデータには、直近のワクチン接種日、抗体価を検査した日、抗体価、抗体検査の方法、などが含まれる。
【0079】
ステップＳ２７において、入力したデータから抗体価予測モデルによる演算を行う。ステップＳ２７では、追加接種またはそれ以降の接種予定者の抗体価の予測分布と予測区間が算出される。予測区間としては、例えば、９５％予測区間を用いることができる。
【0080】
ステップＳ２８において、追加接種またはそれ以降の接種の予測時期を出力する。追加接種またはそれ以降の接種の予測時期としては、予測区間が抗体価の閾値５０を下回る点を追加接種またはそれ以降の接種の予測時期とする。なお、閾値の値は５０に限定されず、適宜選択することが可能である。
【0081】
上述のフローチャートは、例えば、サーバとクライアントデバイスを含むネットワークシステムにおいて実施することが可能である。図１８は、本発明が適用されるネットワークシステムの一例を示す図である。ネットワークシステム２０は、サーバ２１とクライアントデバイス２３とデータベース２２を含む。サーバ２１とクライアントデバイス２３とデータベース２２は、ネットワーク２４を通じて接続している。ネットワーク２４には、任意の通信機器が接続しうる。
【0082】
サーバ２１は、出力装置、入力装置、記憶装置、制御装置、演算装置を備えるコンピュータである。クライアントデバイス２３から受けた依頼と要求（リクエスト）に対して、サーバ２１は応答（レスポンス）をしてサービスを提供する。また、サーバ２１は、データベース２２から検査データの提供を受け、検査データに基づいて抗体価の予測モデルを構築する。図１８に示されたサーバ２１は１つであるが、複数のサーバを備える構成であってもよい。
【0083】
データベース２２は、検査データを蓄積する記憶装置である。ネットワーク２４を通じて、検査機関および医療機関から検査データを受信し、蓄積する。また、データベース２２は、サーバ２１から検査データを提供することもできる。また、データベース２２は、クライアントデバイス２３からワクチンの接種記録や抗体価の検査データを受信し蓄積することもできる。図１８に示されたデータベース２２は１つであるが、複数のデータベースを備える構成であってもよい。
【0084】
クライアントデバイス２３は、スマートフォン、パーソナルコンピュータなどの情報通信機器を含む。任意の数のクライアントデバイスが、ネットワーク２４を通じて接続しうる。
【0085】
ネットワークは、イーサネット、企業規模のコンピュータネットワーク、イントラネット、インターネット、ローカルエリアネットワーク（ＬＡＮ）、広域ネットワーク（ＷＡＮ）、または他のコンピュータネットワークを含むことができる。任意のネットワークを含むことができる。
【0086】
クライアントデバイス２３に追加接種またはそれ以降の接種の予測時期を出力し表示する場合、サーバ・クライアント方式や、クライアントデバイス２３のスタンドアロン方式とすることが考えられる。サーバ・クライアント方式の場合、クライアントデバイス２３がリクエストとしてステップＳ２６を行い、サーバ２１がそのほかのステップＳ２１～Ｓ２５、Ｓ２７をおこなう。スタンドアロン方式とする場合、サーバ２１によって算出された事後分布（ステップＳ２５）をクライアントデバイス２３へ送付し、クライアントデバイス２３がステップＳ２６～ステップＳ２８を行う。
【0087】
図１９は、クライアントデバイス２３の表示例を示す図である。クライアントデバイス２３として、タッチパネルを持つスマートフォンを一例として示す。
【0088】
図１９（ａ）を参照して説明する。入力フィールド２３１は、ワクチン接種日の入力を受け付ける。追加接種またはそれ以降の接種予定者は、ワクチン接種日を入力する。入力フィールド２３２は、抗体価を検査した日の入力を受け付ける。入力フィールド２３３は、抗体検査の方法を入力する。たとえば、検査方法ごとに抗体価の測定単位が異なる場合、抗体検査の方法の情報に基づいて、抗体値の単位が変換される。入力フィールド２３４は、検査で判明した抗体価の入力を受け付ける。入力フィールド２３１から２３４に入力したうえで実行ボタン２３５をタップすると、テキスト表示２３６に追加接種またはそれ以降の接種の予測時期が表示される。
【0089】
【0090】
図１９（ｂ）において、抗体価の推移がグラフ表示２３７に示される。抗体価の推移を視覚的に示し、ユーザにより認識しやすくする。
【0091】
ここではスマートフォンを例示したが、クライアントデバイスはこれに限られない。パーソナルコンピュータやタブレットＰＣでもよい。また、入力フィールド２３１から２３４への入力は、タッチパネルインターフェースに限定されない。接種記録が保存されたＩＣカードや記録媒体から必要な情報を入力することとしてもよい。
【0092】
（作用・効果）
従来、追加接種やそれ以降の接種の時期は年齢や既往の有無、居住地ごとに定められており、必要な時期に追加接種やそれ以降の接種をすることができなかった。また、追加接種やそれ以降の接種の時期になると希望者が殺到し、所望のタイミングで追加接種やそれ以降の接種をすることができなかった。これに対し、本発明の適用例では、個人ごとに高い予測精度で抗体価が閾値を下回る時期を予測し、予測時期を提示することができる。これにより、必要な時期に追加接種やそれ以降の接種を受けることが可能となり、また、追加接種やそれ以降の接種の時期が集中するのを緩和されることも期待される。
従来、追加接種やそれ以降の接種の時期は年齢や既往の有無、居住地ごとに定められており、必要な時期に追加接種やそれ以降の接種をすることができなかった。また、追加接種やそれ以降の接種の時期になると希望者が殺到し、所望のタイミングで追加接種やそれ以降の接種をすることができなかった。これに対し、本発明の適用例では、個人ごとに高い予測精度で抗体価が閾値を下回る時期を予測し、予測時期を提示することができる。これにより、必要な時期に追加接種やそれ以降の接種を受けることが可能となり、また、追加接種やそれ以降の接種の時期が集中するのを緩和されることも期待される。
【0093】
以上、本発明の実施の形態について説明したが、本発明は、上述した実施の形態に限定されるものではなく、本発明の要旨を逸脱しない範囲において種々の変更が可能である。
【符号の説明】
【0094】
２０ ネットワークシステム、２１ サーバ、２２ データベース、
２３ クライアントデバイス、２４ ネットワーク、
２３１～２３４ 入力フィールド、２３５ 実行ボタン、
２３６ テキスト表示、２３７ グラフ表示
２３ クライアントデバイス、２４ ネットワーク、
２３１～２３４ 入力フィールド、２３５ 実行ボタン、
２３６ テキスト表示、２３７ グラフ表示
【図1】
【図2】
【図3】
【図4】
【図5】
【図6】
【図7】
【図8】
【図9】
【図10】
【図11】
【図12】
【図13】
【図14】
【図15】
【図16】
【図17】
【図18】
【図19】