JP2023018569A - 泌乳量算出システム - Google Patents

Abstract

【課題】収集した泌乳量に係る情報が、同一範囲にそろわない場合であっても、要因分析を可能にする泌乳量算出システムを提供する。【解決手段】泌乳量算出システム１は、分娩後日数、実乳量及び気候データを含む要因パラメータを牧場データとして蓄積するメモリ１４と、牧場データの中の一定期間の分娩後日数と実乳量を補間する補間式を作成する補間式作成部と、補間式から算出される特定の分娩後日数の補間乳量を目的変数とし牧場データから選ばれた少なくとも１つの要因パラメータを説明変数として推定泌乳量を算出する回帰式を作成する回帰式作成部と、入力された前記要因パラメータの値に対して回帰式から推定泌乳量を算出する乳量算出部を含む本体１２（制御装置）と、を備える。【選択図】図１

Description

本発明は、分娩後日数に応じた乳牛の搾乳量を、過去のデータに基づいて推定する泌乳量算出システムに関するものである。

酪農業にとっては、泌乳量の増減は、牧場経営に直接影響する大きな問題である。より多くの乳量を得るために、飼料の配合を工夫したり、飼育環境を改善し、乳牛のストレスを少なくするといった工夫が行われている。

一方、このような改善は、個々の牧場に依存する部分が多いため、他の牧場との比較は行いにくい。一般に乳牛は暑さに弱いため、緯度の高い地方の方が泌乳量は多いと考えられるが、気候の違いがどの程度泌乳量に影響しているのかは明らかになっているわけではない。

まして、泌乳量に影響を及ぼす要因は非常に多く考えられるため、どのようなケースの時に、どのような要因が泌乳量に大きく影響を及ぼしているかについては、経験に頼っているのが現状である。

したがって、考えられる要因を入力すると泌乳量を推定できるシステムは、酪農経営に大きな指針を与えるものである。

特許文献１は、情報端末を用いて、乳牛を泌乳ステージ別若しくは繁殖ステージ別に群で把握するシステムを提供している。

また、特許文献２では、牛個体ごとに、過去乳量から回帰式により泌乳曲線を算出するものが開示されている。

特開２０２０－１５６３５９号公報 特開２０２０―０２０７０７号公報

特許文献１は、乳牛の現在の状況を把握しようとするものであり、特許文献２は、過去の泌乳量から未来の泌乳量を推定しようとするものである。しかし、乳牛の泌乳量を増加させるための要因を調べることができるものではない。

高温やストレスによるオキシトシンの減少とアドレナリンの分泌が泌乳量に影響があるという生理的な理屈は理解されているものの、どういった飼育環境や乳牛自身の状態が、乳牛にとって快適な状況であるのかは、人間には理解しにくい。そのため、多くの事象（情報）を集めて要因分析し、泌乳量を高めようとする方法は有用であると考えられる。

ここで、多くの情報とは、できるだけ多くの牧場の、飼育環境や飼料配合、乳牛自身の個別情報が該当する。これらの情報は日々更新されている。また、情報の多様性を確保するには、異なる地域の情報もあれば好ましい。したがって、多くの牧場からのデータを利用できるようにするのが好ましい。

しかしながら、規模や事情の異なる複数の牧場からのデータの定期的な収集は困難を極める。牧場によっては、毎日各種データが数値化される場合もあれば、数日毎でなければ、データがそろわない場合もあるからである。

特に所謂泌乳曲線は分娩後日数と乳量の関係を表すものであるが、分娩後日数が１日違いでそろっている乳牛がいる牧場はよほど大規模でないかぎり、想定しにくい。したがって、実際に収集できる泌乳量のデータでさえ、抜けが生じたり、情報の範囲（分娩後日数の範囲）が異なり、画一的に分析できないという課題があった。

本発明に係る泌乳量算出システムは、上記の課題に鑑みて想到されたものであり、収集した泌乳量に係る情報が、同一範囲にそろわない場合であっても、要因分析を可能にするシステムを提供するものである。

より具体的に本発明に係る泌乳量算出システムは、
分娩後日数、実乳量および気候データを含む要因パラメータを牧場データとして蓄積するメモリと、
前記牧場データの中の一定期間の前記分娩後日数と前記実乳量を補間する補間式を作成する補間式作成部と、
前記補間式から算出される特定の前記分娩後日数の補間乳量を目的変数とし、前記牧場データから選ばれた少なくとも１つの前記要因パラメータを説明変数として推定泌乳量を算出する回帰式を作成する回帰式作成部と、
入力された前記要因パラメータの値に対して前記回帰式から推定泌乳量を算出する乳量算出部を含む制御装置を有することを特徴とする。

本発明に係る泌乳量算出システムは、収集した泌乳量に関するデータを、補間することで、一定の範囲の泌乳量を１日単位でさえ得ることができるので、泌乳量に影響を及ぼす要因の推定に役立つという効果を与える。

また、気温といった予測可能な要因が泌乳量に及ぼす影響が判れば、将来の泌乳量の予測を行うことができるという効果も奏する。

本発明に係る泌乳量算出システムの構成を示す図である。 泌乳量算出システムの全体（メイン）フローを示す図である。 泌乳量算出システムの算出する出力を例示する図である。 推定泌乳量ＥＹを算出する工程のフローを示す図である。 補間式を説明する図である。 補間式を作成するフローを示す図である。 回帰式の作成を説明する図である。 回帰式を作成する工程のフローを示す図である。 図３（ａ）～（ｅ）の推定泌乳量ＥＹに係る曲線を拡大した図である。 図４のステップＳ２０６の要因分析の工程の処理を示す図である。 要因分析において、一定期間ｍの乳量の実測値ＲＹと回帰式から求めた推定泌乳量ＥＹを比較する様子を示す図である。 個別推定泌乳量ＩＹの処理フローを示す図である。 個別推定泌乳量ＩＹを算出する工程の表示例を示す図である。 群推定泌乳量ＧＹを算出する工程の処理を示す図である。 群推定泌乳量ＧＹを算出する工程の表示例を示す図である。 推定泌乳量比較ＣＹの処理フローを示す図である。 推定泌乳量比較ＣＹを算出する工程の出力例を示す図である。 推移泌乳量比較ＴＹを算出する工程の処理フローを示す図である。 推移泌乳量比較ＴＹを算出する工程の出力例を示す図である。 推移泌乳量比較ＴＹを算出する工程の出力例に乳牛個体の実乳量ＲＹを同時にプロットした場合の結果について示す図である。 泌乳量算出システムを実際の牧場に使用し、推定泌乳量は、ＴＨＩと日照時間を要因パラメータとして重回帰できることを示す図である。

以下に本発明に係る泌乳量算出システムについて図面を示し説明を行う。なお、以下の説明は、本発明の一実施形態および一実施例を例示するものであり、本発明が以下の説明に限定されるものではない。以下の説明は本発明の趣旨を逸脱しない範囲で改変することができる。

図１に本発明に係る泌乳量算出システムの構成を示す。本発明に係る泌乳量算出システム１は、端末１０と本体１２およびメモリ１４で構成される。端末１０は、それ自身がＣＰＵ（Ｃｅｎｔｒａｌ Ｐｒｏｃｅｓｓｏｒ Ｕｎｉｔ）を有するコンピュータとメモリおよび表示画面で構成することができ、移動通信端末（所謂スマートフォン）を含めることができる。

本体１２はＣＰＵで構成される。１台のＣＰＵであってもよいし、複数のＣＰＵを接続したものであってもよい。本体１２はサービスを提供する装置なので、サーバーと呼んでよい。メモリ１４は本体１２が使用するメモリである。

本体１２と端末１０は、通信可能に接続されている。通信は、有線であっても、無線であってもよい。端末１０は、クライアントと呼ぶこともできる。端末１０は、実際に乳牛を飼育している牧場１６Ａに設置されるのが好ましい。ここで、「設置」とは、端末１０を牧場１６Ａに物理的に設置する意味の他、牧場１６Ａの責任者若しくは担当者が移動通信端末を保有することであってもよい。もちろん、端末１０は、牧場１６Ａ以外の場所に設置されてもよい。また、牧場１６は、乳牛が飼育されていればよい。

本体１２には、複数の牧場１６が接続されていてもよい。泌乳量算出システム１は、泌乳量に係るデータを多量に集め、これらを分析することで、泌乳量を算出するので、多くのデータが集まるのが望ましいからである。また、本体１２は、外部情報１８とインターネットを介して接続されていてもよい。牧場１６から得られるデータだけでなく、全国的な天気に関する情報をネット経由で取得して利用することもできるからである。

このように、泌乳量算出システム１は、多量のデータを様々な場所から収集し、これらを分析することで、泌乳量を算出するので、ネットワーク経由でサービスを提供する形態が好適な実施形態としてあり得る。すなわち、泌乳量算出システム１はクラウドで構成されてもよい。

したがって、被乳量算出システムを実施するクラウドの形態としては、代表的なＳａａＳ（Ｓｏｆｔｗａｒｅ ａｓ ａ Ｓｅｒｖｉｃｅ）が好適であるが、ＰａａＳ（Ｐｌａｔｆｏｒｍ ａｓ ａ Ｓｅｒｖｉｃｅ）、ＨａａＳ（Ｈａｒｄｗａｒｅ ａｓ ａ Ｓｅｒｖｉｃｅ）若しくはＩａａＳ（Ｉｎｆｒａｓｔｒｕｃｔｕｒｅ ａｓ ａ Ｓｅｒｖｉｃｅ）といった形態であってもよい。

本発明に係る泌乳量算出システム１は、端末１０から要因パラメータＰを入力すると、分娩後日数Ｎ日の乳牛から搾乳できる推定泌乳量^ＮＥＹ（ｋｇ）を得ることができる。ここで要因パラメータＰとは、気候データ、乳牛自体の体重、飼料摂取量といった個別データ、摂取させている飼料の種類といった飼料データ等を含む牧場データの内、泌乳量を最もよく説明することのできる項目である。この要因パラメータＰは、直接測定可能な項目以外に、直接測定可能な項目のデータを加工したデータであってもよい。この要因パラメータＰは、泌乳量算出システム１内で、探すことができる。

気候データは、日の出時刻、日の入り時刻、気温、湿度、日照時間、日射量、風向、風速等である。また、牛舎内温度、牛舎内湿度を含んでもよい。また、ＴＨＩ（Ｔｅｍｐｅｒａｔｕｒｅ Ｈｕｍｉｄｉｔｙ Ｉｎｄｅｘ：温度湿度指数）等、これらの加工データを含んでもよい。なお、ＴＨＩとは一例として（１）式で表される指数である。

また、気候データは、牛舎内、牛舎外の気候だけでなく、牧場の１ｋｍ四方の気候、地域の気候また全国の天気予報サービスからの気候データであってもよい。

個別データは、乳牛自体のデータである。具体的には、登録番号、血統、産回数、分娩後日数、乳量、体重、病歴、所有者、飼育場所等である。

飼料データは、過去に摂取させた飼料の種類や配合、回数などが含まれる。

牧場１６に設置された端末１０は、牧場の乳牛に係る個別データ、実際に摂取させた飼料データ、日々の気候データを牧場データＦＤとして本体１２に送信する。本体１２はこれらのデータをメモリ１４に蓄積し、端末１０から入力された要因パラメータＰに従って、分娩後日数Ｎに応じた推定泌乳量^ＮＥＹを算出し、端末１０に表示させる。

なお、本体１２は、気候データを外部情報１８から取得してもよい。クローズタイプの牛舎であれば、牛舎内の環境は管理されており、データも送信しやすい。しかし、オープンタイプの牛舎の場合、その日の牛舎内の風速は取得しにくい。そのような場合は、外部の気候データを参考にできるからである。

なお、気候に関するデータは、さまざまな物が提供されており、それらも有効に使用することができる。また、牧場データＦＤは端末１０以外の方法で本体１２に送信されてもよい。もちろん、牧場データＦＤを本体１２に送信する端末１０と、要因パラメータＰを入力し、泌乳量を算出させる端末１０が別々であってもよい。

図２に泌乳量算出システム１の全体（メイン）フローを示す。なお、牧場データＦＤを本体１２に送信するのは、各牧場が適宜行うとする。泌乳量算出システム１は、計算条件や要因パラメータＰといった端末等から入力された入力に対して、推定泌乳量ＥＹを算出する工程、個別推定泌乳量ＩＹを算出する工程、群推定泌乳量ＧＹを算出する工程、推定泌乳量比較ＣＹを算出する工程、推移泌乳量比較ＴＹを算出する工程による結果を返すことができる。この結果は端末１０に返してよい。これらを算出する工程は、メインフローのメニューと呼ぶことができる。

なお、本発明は構成の大部分がソフト的に実施される。そこで、処理の「工程」とは、処理のまとまりをいい、本体１２には、「工程」を実行する「部」が存在すると考えてよい。具体的には、本体１２（制御装置）には、推定泌乳量算出部、個別推定泌乳量算出部、群推定泌乳量算出部、推定泌乳量比較算出部、推移泌乳量比較算出部が存在すると言える。

また、これらの処理工程の中に下位の処理のまとまりがある場合は、それも「部」と呼んでよい。例えば、推定泌乳量算出部には、回帰式を作成する工程が１まとまりの処理としてある。これを「回帰式作成部」と呼ぶことができる。さらにこの中には、「補間式を作成する工程」が含まれ補間式作成部と呼んでよい。また、得られた回帰式に要因パラメータＰの値が代入された際に、推定泌乳量を具体的に算出する工程（ステップＳ２１０）等も乳量算出部と呼ぶことができる。また、泌乳量算出システム１は少なくとも推定泌乳量ＥＹを算出する工程（推定泌乳量算出部）を有していればよい。

泌乳量算出システム１がスタートすると（ステップＳ１００）、終了判断を行う（ステップＳ１０２）。終了判断は、端末１０からの終了指示や、本体１２と端末１０の通信の切断であってよい。終了する場合（ステップＳ１０２のＹ分岐）は、終了する（ステップＳ１０４）。終了しない場合（ステップＳ１０２のＮ分岐）は、次に処理を進める。

次に推定泌乳量ＥＹを算出するか、個別推定泌乳量ＩＹを算出するか、群推定泌乳量ＧＹを算出するか、推定泌乳量比較ＣＹを算出するか、推移泌乳量比較ＴＹを算出するかを選択する（ステップＳ１０６、ステップＳ１０８、ステップＳ１１０、ステップＳ１１２、ステップＳ１１４）。それぞれの選択肢では、選択する場合（Ｙ分岐）は、それぞれの処理に移り、選択しない場合（Ｎ分岐）は処理を次に移す。ステップＳ１１４の推移泌乳量比較ＴＹの算出を行わない場合（ステップＳ１１４のＮ分岐）は、終了判断（ステップＳ１０２）に戻る。

各選択ステップ（ステップＳ１０６、ステップＳ１０８、ステップＳ１１０、ステップＳ１１２、ステップＳ１１４）を選択した場合（Ｙ分岐）は、それぞれの工程が行われ、終了判断（ステップＳ１０２）に戻る。具体的には、推定泌乳量ＥＹを求める場合（ステップＳ１０６のＹ分岐）では、推定泌乳量ＥＹを算出する工程が行われる（ステップＳ１１６）。

個別推定泌乳量ＩＹを求める場合（ステップＳ１０８）は、個別推定泌乳量ＩＹを算出する工程が行われる（ステップＳ１１８）。群推定泌乳量ＧＹを求める場合（ステップＳ１１０）は、群推定泌乳量ＧＹを算出する工程が行われる（ステップＳ１２０）。推定泌乳量比較ＣＹを求める場合（ステップＳ１１２）は、推定泌乳量比較ＣＹを算出する工程が行われる（ステップＳ１２２）。また、推移泌乳量比較ＴＹを求める場合（ステップＳ１１４）は、推移泌乳量比較ＴＹを算出する工程が行われる（ステップＳ１２４）。各処理工程を終了したら再び終了判定（ステップＳ１０２）に戻る。

図３には、それぞれの出力例を示す。図３（ａ）は推定泌乳量ＥＹの出力例である。特定の条件（例えば温度）の下で、分娩後日数に対する泌乳量を示すことができる。図では温度がＴ１℃の時と、Ｔ２℃の時の推定泌乳曲線Ｍが示されている。なお、本発明の泌乳量算出システム１が提供する推定泌乳曲線Ｍは、最小で１日単位の算出値（推定泌乳量ＥＹ）の集合である。これらの算出値は直線もしくは曲線で結ばれていてもよい。推定泌乳量ＥＹで形成した推定泌乳曲線を符号Ｍで表す。言い換えると、推定泌乳曲線Ｍは推定泌乳量ＥＹを適当な分娩後日数Ｎ毎にプロットし、直線若しくは曲線で繋げたものである。

図３（ｂ）は個別推定泌乳量ＩＹの出力例である。個別推定泌乳量ＩＹは、指定された母集団の推定泌乳曲線Ｍに個別の乳牛の実乳量ＲＹをプロットしたものである。図３（ｂ）では、黒丸ｄで表した。実乳量ＲＹが推定泌乳量ＥＹに対して、一定以上多い個体や、一定以上少ない個体は、特にわかるように表示することができる（図３（ｂ）では、符号ｄｕ、ｄｄで示した。）。また、このような個体があった場合は、所有者に通知を行ってもよい。

図３（ｃ）は群推定泌乳量ＧＹの出力例である。群推定泌乳量ＧＹは、指定された母集団の推定泌乳曲線Ｍに群で管理しているグループの乳牛の実乳量ＲＹをプロットしたものである。同じクループの乳牛は同じ印でプロットする。図３（ｃ）では、ユニットＵ１とユニットＵ２の２つの群に属する乳牛に対して、ユニットＵ１は、白丸。ユニットＵ２は黒丸で表した。

また、実乳量ＲＹが推定泌乳量ＥＹに対して、一定以上高い、若しくは一定以上低い場合は、個別推定泌乳量ＩＹの場合同様に、通知を行ってもよい。なお、同一グループであっても、推定泌乳量ＥＹから一定以上離れる乳量に対しては、別のマークを用いてもよい。

図３（ｄ）は推定泌乳量比較ＣＹの出力例である。推定泌乳量比較ＣＹは、例えば環境の異なる牧場の牛の推定泌乳曲線Ｍ１とＭ２等を比較する。なお、この推定泌乳曲線に特定の乳牛の実乳量ＲＹを追加して表示してもよい。

図３（ｅ）は、推移泌乳量比較ＴＹの出力例である。ここでは、同一牧場の年毎８月の推定泌乳曲線を比較したものである。もちろん、比較の対象は、他の牧場との比較であってもよい。また、後述するように、この推定泌乳曲線Ｍに特定の乳牛の実乳量ＲＹを追加してプロットしてもよい。その乳牛の時間的な実乳量の変移が、推定泌乳量ＥＹに対して、どういう位置づけになるかを示すことができる。

図２を再度参照する。以上の全体フローは、本体１２と端末１０の間で分担して行ってよい。例えば、各処理の選択および最終結果の表示は端末１０で行い、各処理は本体１２で行うなどである。以下、各工程を詳説する。

＜推定泌乳量ＥＹ＞
図４に推定泌乳量ＥＹの算出工程のフローを示す。推定泌乳量ＥＹの算出工程が始まったら（ステップＳ１１６）、母集団データの選択と作成条件を決める（ステップＳ２００）。図１で示したメモリ１４中には、多数の牧場データが蓄積されている。その中で、必要な場所およびデータの期間を選択する。また、推定泌乳量ＥＹを算出する条件もここで入力する。

例えば、場所は、自分の牧場で、データ期間は昨年度１年等である。また、条件とは、推定泌乳量ＥＹを算出するための条件すべてを含む。例えば、牧場データを参照する際にまとめる期間の長さを決める等である。より具体的にいうと、実乳量は１日毎の乳量にするか、３日毎、１週間毎、１か月毎の平均乳量であるとか、１か月の最大乳量値である等である。また、これに伴い、その他の変数も同様に一定期間毎に代表値を決めてよい。

次に補間式を作成する（ステップＳ２０２）。図５に補間式について説明する。図５は横軸が分娩後日数Ｎ（日）であり、縦軸は実乳量ＲＹ（ｋｇ）である。例えば、ある牧場の複数の乳牛の、ある期間（例えば２０２０年８月１日）のデータである。ここでプロットされた各点は、乳牛毎の実乳量ＲＹなので、散布図となる。しかし、実乳量であるために、プロットできない範囲も生じる。該当する分娩後日数の乳牛がいない場合である。図５では、無データ領域ＶＲで表した。

このように、牧場データに抜けがあると後の回帰式を作成する際のデータに不足が生じ、回帰の精度が低下する。そこで、この散布図を適当な関数で近似する。好適に利用できるのはＷＯＯＤ曲線である。ＷＯＯＤ曲線は、（２）式で表される曲線で、分娩後日数Ｎに対する泌乳量Ｙを表す曲線としてよく知られている。なお、実際の散布図を補間するには、ＷＯＯＤ曲線に限定される必要はなく他の関数であってもよい。

なお、ここでＡ、Ｂ、Ｃは定数であり、Ｙは泌乳量、Ｎは分娩後日数である。このＷＯＯＤ曲線を最小二乗法で散布図にフィットするように定数Ａ、Ｂ、Ｃを決めることができる。このように、分娩後日数Ｎと乳量Ｙの関係を連続関数でフィッテングさせることを「補間式を作成する」という。散布図を補間した式をより一般化して（３）式と表す。

つまり（３）式の補間式は、乳量Ｙは分娩後日数Ｎの関数として表される。関数の形としては、ＷＯＯＤ曲線は好適に利用できるが、これに限定されるものではない。図４の補間式作成工程（ステップＳ２０２）は、上記の様に実乳量ＲＹの散布図から補間式（３）を決める工程である。

［補間式作成］
図６に補間式を作成する工程（ステップＳ２０２）の処理を示す。補間式を作成する工程が始まると（ステップＳ２０２）、図４（ＥＹ算出フロー）のステップＳ２００で決めた母集団と条件より分娩後日数Ｎと実乳量ＲＹのデータを取り出す（ステップＳ２３０）。これらは実データと呼ぶことができる。具体的に散布図を描いてもよい。実データとは、各乳牛毎の分娩後日数Ｎと実乳量ＲＹの組である。

次に分娩後日数Ｎと乳量Ｙの実データを最もよく反映する補間式を求める。例えば、上記に示した、ＷＯＯＤ曲線を使った最小二乗法を適用して、定数Ａ、Ｂ、Ｃを決定する（ステップＳ２３２）（図５も参照）工程である。これは、実乳量ＲＹにフィッテングさせる関数を求めると言える。そして、（３）式の補間式を得る（ステップＳ２３４）。その後ＥＹ算出ルーチン（図４）に戻る（ステップＳ２３６）。

補間式は、毎日分作成することができる。搾乳はほぼ毎日行われるからである。しかし、環境や乳牛の各個体について、大きな変化がない場合は、３日毎や、１週間毎の実データの平均を実データとみなして補間式（３）を作成してもよい。

また、時間軸を大きくとらえ、１か月の実データの平均値を実データとみなして補間式を作成してもよい。より具体的には、ある乳牛個体の１か月の乳量の平均をその乳牛の実乳量ＲＹとする。またその乳牛の分娩後日数はその１か月の平均（つまり、月初めにＮ＝１０日であれば、その月の分娩後日数を２５日とする）とする等である。これらは、図４の「母集団データ、作成条件」（ステップS２００）で決められていても良い。

引き続き図４を参照する。補間式を作成したら要因分析を行うか否か判断する（ステップＳ２０４）。この判断は泌乳量算出システム１の利用者が端末１０から指示することができる。要因分析とは、牧場データの中から（２）式の補間式から得られる乳量Ｙを最もよく説明できる要因パラメータＰを調べるか否かの判断である。例えば、泌乳量算出システム１を初めて使用する場合や、母集団を大きく変える場合等に利用する。

要因分析を行う場合（ステップＳ２０４のＹ分岐）は、要因分析の工程を行う（ステップＳ２０６）。要因分析工程の詳細は後述する。行わない場合（ステップＳ２０８のＮ分岐）は、処理を次に移す。要因分析については、図１０で詳説する。

次に回帰式を作成する（ステップＳ２０８）。図７に回帰式の作成について説明する。図５、図６で示した補間式（３）（図７（ａ））によって、例えば各月毎の分娩後日数Ｎの乳牛の補間乳量が算出できる。ここでは補間式を一般形Ｆｗ（）で表し、乳量Ｙは月ごとの平均乳量とする。したがって、Ｙ_１月は１月の補間乳量を表し、分娩後日数Ｎの乳牛の補間式（３）を求める式をＦｗ_１月（Ｎ）と表す。

ここで、分娩後日数が５０日後をＮ_５０と表す。１月から１２月までの各月における分娩後日数が５０日後の乳牛の補間乳量は、Ｆｗ_１月（Ｎ_５０）、Ｆｗ_２月（Ｎ_５０）、・・・、Ｆｗ_１２月（Ｎ_５０）で算出される。もちろんＦｗ_ｍ（Ｎ）は、各月毎に補間式が作成されている（「ｍ」は「月」を表す。）。

この各月の補間乳量を最もよく説明することのできる因子を要因パラメータＰとして決める。どのような要因パラメータＰがよいのかは、ステップＳ２０６の要因分析の工程を行って知ることになる。すでに推定泌乳量ＥＹを求めるために使用する要因パラメータＰが決まっている場合は、それを使用する。要因パラメータＰは１つでなくてもよい。すなわち、これは各月ごとの補間乳量Ｙを１以上の因子によって最小二乗法による回帰式を求めることである。

よく知られているように回帰式は（４）式のように表される。

ここで、ＥＹは推定泌乳量であり、ｘ_１、ｘ_２、・・・、ｘ_ｋは因子（要因パラメータＰ）であり、ａ_１、ａ_２、・・・、ａ_ｋ、ｃは定数である。

図７（ｃ）には、因子が１つの場合の回帰式の例を示す。図７（ａ）で示す補間式（３）から、分娩後日数５０日の乳牛の各月毎の補間乳量Ｙ_ｍを算出する（図７（ｂ）参照。）。ここで、ｍは月を表す。各月毎の補間乳量とは、月初めの３日間の平均とするなどであってもよい。特定の１日だけのデータにその月の補間乳量を代表させてもよい。

次に、これら月ごとの補間乳量を要因パラメータＰ（ここでは、例えば正午の温度）で整理すると図７（ｃ）のグラフを得ることができる。なおここでは、各月の分娩後５０日の乳牛の推定泌乳量ＥＹは、正午の気温でよく説明できるとした場合の結果を例示する。回帰式は、（５）式で表される。

ここで、ｘ_１は、正午の温度であり、ＥＹは推定泌乳量である。なお、このような回帰式は分娩後日数分だけ作成される。つまり、分娩後日数Ｎの最終日を３００日とすると、１日から３００日後までの回帰式を得ることができる（図７（ｅ）参照。）。

因子が複数ある場合は、図７（ｃ）のように２次元では記載できなくなるが、複数の要因パラメータによる重回帰式を得れば、分娩後日数Ｎと補間乳量Ｙの関係をよく近似することができる場合がある。このようにして求めた回帰式を（６）式と表す（図７（ｄ）参照。）。

ここで、^ＮＥＹは、分娩後日数Ｎの乳牛が要因パラメータＰの時の推定泌乳量ＥＹと呼ぶ。

［回帰式作成］
図８に回帰式を作成する工程のフローを示す。図４で示したＥＹ算出フローで回帰式を作成する工程（ステップＳ２０８）が開始されると、図８に飛び、まず要因パラメータＰが入力される（ステップＳ２５０）。要因パラメータＰは、端末１０から入力することができる。つまり、泌乳量算出システム１の利用者が入力する。ここでは仮に要因パラメータＰがｐ_１（温度）であったとする。

次に要因パラメータｐ_１について指定された一定期間ｍ（これはステップＳ２００で定められる）毎に、分娩後日数Ｎ毎の推定泌乳量Ｙｍを、補間式（３）によって求める（ステップＳ２５２）。つまり、図７（ｂ）の補間乳量Ｙｍを求める。これによって、図７（ｃ）の散布図を描くことができる。

次に要因パラメータｐ_１を説明変数、補間乳量Ｙｍを目標変数として、回帰式を求める（ステップＳ２５４）。すなわち、この時要因パラメータｐ_１について、分娩後日数Ｎが同じ場合の乳量の関係（６）式（図７（ｄ））を得ることができる。回帰式の求め方は公知の方法でよい。

このように本発明では、要因パラメータＰを使って推定泌乳量を算出する際に、実乳量ＲＹを直接使用するのではなく、実乳量ＲＹを補完する補間式から算出した補間乳量を利用する。そのため、実乳量ＲＹにデータの抜けがあっても、妥当な値を平均的な乳量として算出することができる。

（６）式は、分娩後日数Ｎが１から最終日（ここで最終日は、分娩後日数の最長日の意味）までの数を作ることができる（図７（ｅ）参照。）。回帰式（６）は、図７（ｅ）をまとめて表している。回帰式が求まったら、図４の推定泌乳量ＥＹを求めるルーチンに戻る（ステップＳ２５６）。

図４を再び参照し、推定泌乳量ＥＹを求める回帰式（６）が求められれば、回帰式（６）に要因パラメータＰを入力すれば推定泌乳量ＥＹを求めることができる（ステップＳ２１０）。

分娩後日数Ｎに対応する回帰式（６）は要因パラメータＰの値を入力すると、推定泌乳量ＥＹを算出する。したがって、分娩後日数Ｎ毎を表す回帰式（６）に順次同一の要因パラメータＰの値を入力することで、同一の要因パラメータＰの値における分娩後日数Ｎと推定泌乳量ＥＹのデータセットを得ることができる。図３（ａ）～（ｅ）の推定泌乳量に係る曲線は、このようにして得られた分娩後日数Ｎと推定乳量ＥＹのデータセットをプロットし、若しくはこの点同士を直線でつないだものである。

図９に図３（ａ）～（ｅ）の推定泌乳量ＥＹをプロットした推定泌乳曲線Ｍを拡大した図を示す。分娩後日数毎に推定泌乳量ＥＹをプロットしたものである。要因パラメータＰは、ｐ_１とｐ_２の２つの要因で重回帰式を求めた場合を示している。分娩後日数Ｎは、飛び飛びの値（１０日、２０日、３０日、５０日、８０日、１００日、１５０日、２００日、２５０日）で求めたが、１日単位で求めることもできる。各推定泌乳量ＥＹは、（６）式が表す別々の回帰式（分娩後日数違い。図７（ｅ）参照。）から求められたものである。すなわち、これらの式の定数は異なる。

推定泌乳量ＥＹを求めるこのステップＳ２１０（図４参照）は、クライアント側が端末１０から入力し、本体１２が推定泌乳量ＥＹを返す。クライアント側は、少なくとも１つの推定泌乳量ＥＹを返せば本発明に係る泌乳量算出システム１としての要件を満たす。複数の分娩後日数Ｎに対して推定泌乳量ＥＹを求め、図９のようにグラフにしてもよい。

図４を再度参照し、要因パラメータＰの再入力を行う場合（ステップＳ２１２のＹ分岐）は、再びステップＳ２１０を行う。推定泌乳量ＥＹを求めない場合（ステップＳ２１２のＮ分岐）は、終了判定を行う（ステップＳ２１４）。終了する場合（ステップＳ２１４のＹ分岐）は、メインルーチンへ戻る（ステップＳ２１６）。終了しない場合（ステップ２１４のＮ分岐）は、ステップ２００に戻り、ＥＹ算出工程を再実行する。

［要因分析］
図１０には、図４のステップＳ２０６の要因分析の工程の処理を示す。要因分析の工程は、回帰式の作成（ステップＳ２０８）と途中までは同じである。具体的には、図８で示すステップＳ２５０、ステップＳ２５２、ステップＳ２５４は、図８の場合と同じである。したがって、ステップ番号も同じ番号を用いる。

したがって、要因パラメータＰを入力し（ステップＳ２５０）、要因パラメータＰについて指定された一定期間ｍ（これはステップＳ２００で定められる）毎に、分娩後日数Ｎ毎の乳量Ｙｍを、補間式（３）によって求め（ステップＳ２５２）、要因パラメータＰを説明変数、乳量Ｙｍを目標変数として、回帰式を求る（ステップＳ２５４）工程が行われる。

要因分析では、一定期間ｍの乳量の実測値ＲＹと回帰式から求めた推定泌乳量ＥＹを比較する（ステップＳ２６０）。図１１はこの様子を図にしたものである。要因パラメータとしては、ｐ_１とｐ_２の２つを選んだ場合を示す。例えば、温度と風速等である。一定期間は異なる４日の場合を示した。例えば、春夏秋冬の代表的な1日としてよい。具体的には、ｍ１月ｎ１日、ｍ２月ｎ２日、ｍ３月ｎ３日およびｍ４月ｎ４日とする。これらの日の具体的な要因パラメータｐ_１、ｐ_２は記録として残っている。要因パラメータの選択は、別途主要因分析を行ってもよいし、試行錯誤によって選択してもよい。要因パラメータｐ_１、ｐ_２を入力すれば回帰式（６）から推定泌乳量ＥＹを算出することができる。

分娩後日数Ｎは１から最大日まで求めることができるが、ここでは、３０日、５０日、１００日、１５０日の４種について求めた状態を示している。これらの値をそれぞれの日時の実乳量ＲＹとの差の絶対値ＥＲＲとする。なお、図１１で^３０ＲＹは、分娩後日数が３０日の乳牛の実乳量である。また、^３０ＥＲＲは、分娩後日数が３０日の場合の推定泌乳量ＥＹと実乳量ＲＹの差の絶対値を示す。このＥＲＲは、分娩後日数が５０日、１００日、１５０日の場合も算出している。なお、後述する実施例では、分娩後日数が５０日の場合について、推定泌乳量ＥＹと実乳量ＲＹの比較を行った例を示した。

４日分について、推定泌乳量ＥＹと実乳量ＲＹの差の絶対値ＥＲＲを合計したものを総誤差ＴΣで表す。このようにステップ２６０では、推定泌乳量ＥＹと実乳量ＲＹを比較する。

再度図１０を参照して、推定泌乳量ＥＹと実乳量ＲＹを比較した後は、終了判定を行う（ステップＳ２６２）。終了判定は、ステップＳ２５０で入力した要因パラメータＰが推定泌乳量ＥＹと実乳量ＲＹの差の絶対値ＥＲＲが妥当であるか否かで決める。ＥＲＲが妥当であるか否かについては、予め値を決めておいてもよいし、利用者が判断してもよい。また、最適値を選択するようにプログラム（所謂「ＡＩ」を含んでよい。）を設定してもよい。また、推定泌乳量ＥＹと実乳量ＲＹの差は絶対値として説明したが、２乗値を用いてもよい。

終了する場合（ステップＳ２６２のＹ分岐）は、ＥＹ算出フローに戻る（ステップＳ２６４）。継続する場合は、再びステップＳ２５０に戻り、要因パラメータＰを入力しなおす。

再度図４を参照する。要因分析（ステップＳ２０６）を抜けたら、乳量算出の工程（ステップＳ２１０）に処理を移す。要因分析の工程（ステップＳ２０６）を抜けたら、好適な要因パラメータＰを見つけていると考えられるので、推定泌乳量ＥＹを算出する工程（ステップ２１０）を行ってよいと考えられるからである。これ以降の工程は上記に説明した通りである。

なお、本発明に係る泌乳量算出システム１においては、回帰式（６）を求める工程は、いずれのメニューでも利用される。したがって、ステップＳ２００、ステップＳ２０２、ステップＳ２０４、ステップＳ２０６、ステップＳ２０８までをまとめてステップＳ２９０とする。ステップＳ２９０は、回帰式取得工程と呼んでもよい。

推定泌乳量ＥＹの算出（ステップＳ１１６）は、要因パラメータＰの値を入力すると、推定泌乳量ＥＹを算出する。この場合の要因パラメータＰの値は、過去のデータである必要はなく、将来の予想する要因パラメータＰの値であってもよい。すなわち、明日の要因パラメータＰの値若しくは半年先の要因パラメータＰの値であってもよい。

以上のように推定泌乳量算出は、特定の母集団に対して、要因パラメータＰが、ある値の時の平均的な乳量を推定した値を返す。言い返すと、推定泌乳量値とは、要因パラメータＰが、ある値の時のその母集団に属する乳牛から見込める予定乳量である。

＜個別推定泌乳量ＩＹ＞
図２を再度参照する。推定泌乳量ＥＹを選択しなかった場合（ステップＳ１０２のＮ分岐）は、個別推定泌乳量ＩＹを算出するか否かが判断される（ステップＳ１０６）。この判断はクライアント側から端末を通じて選択することができる。個別推定泌乳量ＩＹを算出する工程では、推定泌乳量ＥＹと実乳量ＲＹを比較し、その差を算出する。個別推定泌乳量を算出する場合（ステップＳ１０６のＹ分岐）は、個別推定泌乳量ＩＹの算出工程（ステップＳ１１８）に処理が移される。

図１２に個別推定泌乳量ＩＹの算出工程のフローを示す。個別推定泌乳量ＩＹを算出する工程が開始されると（ステップＳ１１６）、回帰式取得工程が行われる（ステップＳ２９０）。これは図４のステップＳ２９０と同じである。したがって、ステップＳ２９０を実行すると、推定泌乳量ＥＹを求める回帰式（６）が得られている。

次に特定の乳牛Ｉの実乳量ＩＲＹを指定する（ステップＳ３００）。この指定は、クライアント側が端末を通じて入力することで行うことができる。実乳量ＩＲＹは、特定の乳牛を指定することで行ってもよい。特定の乳牛の所定期間の実乳量はサーバー側のメモリ１４に記録されているからである。実乳量ＩＲＹを指定すると、その実乳量ＩＲＹを泌乳した乳牛が特定され、その実乳量ＩＲＹを泌乳した時の分娩後日数Ｎおよび要因パラメータＰの値も得ることができる。この時の要因パラメータＰをｐ_１とする。

次に推定泌乳量ＥＹを算出し表示する（ステップＳ３０２）。ここでは、所定の期間の母集団から得られた回帰式（６）から複数の分娩後日数Ｎの回帰式を選択し、要因パラメータＰを入力し、推定泌乳量ＥＹを算出する。分娩後日数Ｎを複数入力することで、分娩後日数と推定泌乳量ＥＹは、２次元表示することができる。入力する分娩後日数Ｎは予め決めておいてもよい。

図１３に個別推定泌乳量ＩＹの表示例を示す。図１３を参照して、横軸は分娩後日数（日）であり、縦軸は乳量Ｙである。要因パラメータＰはｐ_１の場合である。推定泌乳量ＥＹは滑らかな線で表したが、実際は入力した分娩後日数Ｎ個の点をつないだものである。すなわち、推定泌乳曲線Ｍである。

再度図１２を参照する。次に指定された分娩後日数Ｎの乳牛の実乳量^ＮＩＲＹを分娩後日数Ｎの推定泌乳量ＥＹと比較する（ステップＳ３０４）。この比較は、推定泌乳量ＥＹの２次元表示に実乳量ＩＲＹをプロットｄとして表示させてもよい（図１３参照。）。

図１３では、分娩後日数Ｎがｎ_１、ｎ_２、ｎ_３、ｎ_４の４つの点をプロットした。これらのデータは、すべて要因パラメータＰがｐ_１である。また、このように、個別推定泌乳量ＩＹを算出する工程では、複数の個体の実乳量ＲＹをプロットしてもよい。

図１２を再度参照し、実乳量^ＮＩＲＹと推定泌乳量^ＮＥＹの差が優良閾値Ｔｈｕ以上であれば（ステップＳ３０４のＹ分岐）、優良通知を行う（ステップＳ３１０）。優良通知を取得された特定の乳牛の実乳量ＩＲＹは、推定泌乳量ＥＹより多く、優秀な泌乳量（優秀な乳牛）であったことを表している。優良通知は、テキストだけで通知してもよいし、推定泌乳量ＥＹのグラフ上に優良点（ｄｕ）として表示してもよい（図１３参照。）。図１３では、分娩後日数Ｎがｎ２の場合の実乳量ｄｕは、推定泌乳量^ｎ２ＥＹより優良閾値Ｔｈｕ以上大きいことを示している。

実乳量^ＮＩＲと推定泌乳量^ＮＥＹの差が優良閾値Ｔｈｕより大きくなければ（ステップＳ３０４のＮ分岐）、推定泌乳量^ＮＥＹと実乳量^ＮＩＲＹとの比較を行う（ステップＳ３０６）。この比較は、推定泌乳量ＥＹの２次元表示に実乳量ＩＲＹをプロットｄとして表示させてもよい（図１３参照。）。

推定泌乳量^ＮＥＹと実乳量^ＮＩＲＹの差が劣勢閾値Ｔｈｄより大きければ（ステップＳ３０６のＹ分布）、劣勢通知を行う（ステップＳ３１２）。特定された乳牛の泌乳量ＩＲＹは、推定泌乳量ＥＹより少なく、劣勢な泌乳量（劣勢な乳牛）であったことを表している。劣勢通知は、テキストだけで通知してもよいし、推定泌乳量ＥＹのグラフ上に劣勢点（ｄｄ）として表示してもよい（図３（ｂ）参照。）。

図１３では、分娩後日数Ｎがｎ１の場合の実乳量ｄｄは、推定泌乳量^ｎ１ＥＹに対して劣勢閾値Ｔｈｄより小さいこと（推定泌乳量^ｎ１ＥＹは実乳量ｄｄより劣勢閾値Ｔｈｄ以上大きい。）を示している。

なお、推定泌乳量ＥＹとの差が優勢閾値Ｔｈｕ若しくは劣勢閾値Ｔｈｄほど離れていない場合は、通知は行わなくてよい。図１３では、分娩後日数ｎ３およびｎ４の時の場合である。しかし、この場合でも、実乳量^ＮＩＲＹと推定泌乳量^ＮＥＹの差は求められている。

推定泌乳量^ＮＥＹと実乳量^ＮＩＲＹの差が劣勢閾値Ｔｈｄより大きくなければ（ステップＳ３０６のＮ分岐）、実乳量ＩＲＹの再指定か否かを判断する（ステップＳ３０８）。実乳量ＩＲＹを再指定する場合（ステップＳ３０８のＹ分岐）は、ステップＳ３０２に処理を移す。

実乳量ＩＲＹを再指定しない場合（ステップＳ３０８のＮ分岐）は、終了判定に処理を移す（ステップＳ３１４）。終了しない場合（ステップＳ３１６のＮ分岐）は、回帰式取得工程（ステップＳ２９０）に戻る。終了する場合（ステップＳ３１４のＹ分岐）は、メインルーチンに戻る（ステップＳ３１６）。

以上のように個別推定泌乳量ＩＹを算出する工程では、母集団から得られる推定泌乳量ＥＹに対して、指定された乳牛の泌乳量ＩＲＹを比較し、推定泌乳量ＥＹとの差を求めた上で、推定泌乳量ＥＹから一定以上離れた値の泌乳量の乳牛に対しては、通知が行われる。

＜群推定泌乳量ＧＹ＞
図２を再度参照する。個別推定泌乳量ＩＹを選択しなかった場合（ステップＳ１０８のＮ分岐）は、群推定泌乳量ＧＹを算出するか否かが判断される（ステップＳ１１０）。この判断はクライアント側から端末を通じて選択することができる。群推定泌乳量ＧＹを算出する工程では、同一要因パラメータＰの時に、グループ同士の泌乳量を比較する（差が求められる。）。

群推定泌乳量ＧＹを算出する場合（ステップＳ１１０のＹ分岐）は、群推定泌乳量ＧＹの算出工程（ステップＳ１２０）に処理が移される。

図１４に群推定泌乳量ＧＹの算出工程の処理を示す。群推定泌乳量ＧＹを算出する工程が開始されると（ステップＳ１２０）、回帰式取得工程が行われる（ステップＳ２９０）。これは図４のステップＳ２９０と同じである。したがって、ステップＳ２９０を実行すると、推定泌乳量ＥＹを求める回帰式（６）が得られている。この際、母集団には、比較されるグループの乳牛は当然含まれている。

次に群が指定される（ステップＳ３５０）。「群」とは、２頭以上の乳牛のグループで、飼料、環境、産次数等をそろえた乳牛のまとまりである。以後Ａ群とＢ群があり、それぞれ分娩後日数が異なる乳牛が３頭ずつ属しているとする。

群が指定されると、推定泌乳量ＥＹが算出される。もちろん、推定泌乳曲線Ｍを表示させてもよい（ステップＳ３５２）。母集団（ここでは少なくともＡ群とＢ群の乳牛）が決まるので、要因パラメータＰによる推定泌乳量ＥＹが算出できる。

次にＡ群とＢ群に属する乳牛に対して、１頭ずつ実乳量^ＮＩＲＹと推定泌乳量^ＮＥＹが比較される（ステップＳ３５４およびステップＳ３５６）。そして、優秀閾値Ｔｈｕより実乳量ＩＲＹが高ければ優秀通知を行い（ステップＳ３５４）、劣勢閾値Ｔｈｄより実乳量ＩＲＹが低ければ劣勢通知を行う（ステップＳ３６２）。

次にＡ群、Ｂ群の全ての乳牛に対して評価が行われたかを判断し（ステップＳ３５８）、終了していなければ（ステップＳ３５８のＮ分岐）、再度実乳量^ＮＩＲＹと推定泌乳量^ＮＥＹが比較される（ステップＳ３５４およびステップＳ３５６）。

全ての乳牛に対する評価が終了したら（ステップＳ３５８のＹ分岐）、表示を行う（ステップＳ３６４）。図１５に群推定泌乳量ＧＹの表示例を示す。ここでは白丸をＡ群とし、黒丸をＢ群とした。Ａ群は推定泌乳量ＥＹと比較して高めの泌乳量を示し、Ｂ群は推定泌乳量ＥＹと比較して低めの泌乳量を示している。これらのグループは、要因パラメータＰ以外の要因で差がついていると考えられる。なお、優勢通知を受けたものは白抜き星印で表し、劣勢通知を受けたものは黒四角で表した。

以上のように、群推定泌乳量ＧＹでは、指定した群の乳牛の泌乳量が推定泌乳量ＥＹを基準として比較することができる。

＜推定泌乳量比較ＣＹ＞
図２を再度参照する。群推定泌乳量ＧＹを選択しなかった場合（ステップＳ１１０のＮ分岐）は、推定泌乳量比較ＣＹを算出するか否かが判断される（ステップＳ１１２）。この判断はクライアント側から端末を通じて選択することができる。推定泌乳量比較ＣＹとは、推定泌乳量ＥＹ同士を比較する。すなわち、要因パラメータＰの値が違う推定泌乳量ＥＹ同士を比較する。若しくは、母集団同士が異なる場合の比較であってもよい。

推定泌乳量比較ＣＹを算出する場合（ステップＳ１１２のＹ分岐）は、推定泌乳量比較ＣＹの算出工程（ステップＳ１２２）に処理が移される。

図１６に推定泌乳量比較ＣＹの算出工程の処理フローを示す。推定泌乳量比較ＣＹを算出する工程が開始されると（ステップＳ１２２）、回帰式取得工程が行われる（ステップＳ２９０）。ただし、ここでは、母集団違いの回帰式が取得されてもよい。例えば、異なる地方の牧場同士の比較である。これは図４のステップＳ２９０で、母集団を複数選択するのと同じである。もちろん、母集団は１つであってもよい。したがって、ステップＳ２９０を実行すると、推定泌乳量ＥＹを求める回帰式（６）が得られている。ここでは、母集団は符号Ｕとして以下説明する。母集団Ｕは３つ以上あってもよい。

回帰式が取得されると、母集団Ｕおよび要因パラメータＰが指定される（ステップＳ３８０）。母集団Ｕおよび要因パラメータＰが指定されると、推定泌乳量ＥＹが算出され、表示される。

その後再指定の要否が問い合わされ（ステップＳ３８４）、再度母集団Ｕおよび要因パラメータＰを指定する場合（ステップＳ３８４のＹ分岐）は、母集団Ｕおよび要因パラメータＰが指定される処理（ステップＳ３８０）に処理は移され、再指定が不要であれば（ステップＳ３８４のＮ分岐）、メインフローに戻る（ステップＳ３８６）。

図１７には、推定泌乳量比較ＣＹ工程の出力例を示す。図１７（ａ）は、同一母集団Ｕの場合であって、要因パラメータＰが異なる場合である。例えば、同一牧場の季節違いの表示等である。

図１７（ｂ）は、要因パラメータＰが同じであって、母集団がＵ１およびＵ２と異なる場合の結果である。立地条件の違いで、気候の全く異なる牧場同士等で利用することができる。

このように推定泌乳量比較ＣＹ工程では、母集団Ｕ同士が異なる場合の推定泌乳量ＥＹの比較や、同一母集団であっても、要因パラメータＰが異なる場合などの推定泌乳量を算出し表示させることができる。

＜推移泌乳量比較ＴＹ＞
図２を再度参照する。推定泌乳量比較ＣＹを選択しなかった場合（ステップＳ１１２のＮ分岐）は、推移泌乳量比較ＴＹを算出するか否かが判断される（ステップＳ１１４）。この判断はクライアント側から端末を通じて選択することができる。

推移泌乳量比較ＴＹを算出する工程では、時間軸の異なる推定泌乳量ＥＹ同士の比較が行われる。要因パラメータＰは気候データを含むので、異なる要因パラメータＰ同士の推定泌乳量ＥＹの比較も時間軸の異なる推定泌乳量ＥＹ同士の比較と言える。したがって、要因パラメータＰの値が違う推定泌乳量ＥＹ同士を比較するという意味で、推移泌乳量比較ＴＹは、推定泌乳量比較ＣＹと似ている。しかし、推移泌乳量比較ＴＹでは、比較する日時が先に決定され、その日時の要因パラメータＰが使用されるという点で異なる。推定泌乳量比較ＣＹを算出する工程では、要件パラメータＰが先に決定される。

また、推移泌乳量比較ＴＹを算出する工程では、特定の乳牛の実乳量の継時変化をその時の要因パラメータＰとの比較で表すこともできる。

推移泌乳量比較ＴＹを算出する場合（ステップＳ１１４のＹ分岐）は、推移泌乳量比較ＴＹの算出工程（ステップＳ１２４）に処理が移される。

図１８に推移泌乳量比較ＴＹの算出工程の処理フローを示す。推移泌乳量比較ＴＹを算出する工程が開始されると（ステップＳ１２４）、回帰式取得工程が行われる（ステップＳ２９０）。ただし、ここでは、母集団違いの回帰式が取得されてもよい。例えば、同一牧場における異なる年度のデータ同士を比較する場合である。これは図４のステップＳ２９０で、母集団を複数選択するのと同じである。もちろん、母集団は１つであってもよい。したがって、ステップＳ２９０を実行すると、推定泌乳量ＥＹを求める回帰式（６）が得られている。ここでは、年度は符号Ｑで表す。年度Ｑは複数であってもよく、また、数か月毎の期間であってもよい。

回帰式が取得されると、表示期間が指定される（ステップＳ４００）。表示期間とは、例えば各年で表示比較したい月や週、日時であってもよい。決められた母集団において、日時が決まると、要因パラメータＰを決めることができる。また、表示期間には、乳牛自体の指定を含めてもよい。推定泌乳量の変移と共に、個体の泌乳量ＲＹの変移を知りたい場合もあるからである。要因パラメータＰが決まれば、推定泌乳量ＥＹが算出され、表示される（ステップＳ４０２）。

その後再指定の要否が問い合わされ（ステップＳ４０４）、表示期間を再指定する場合（ステップＳ４０４のＹ分岐）は、表示期間が指定される処理（ステップＳ４００）に処理は移され、再指定が不要であれば（ステップＳ４０４のＮ分岐）、メインフローに戻る（ステップＳ４０６）。

図１９には、推移泌乳量比較ＴＹ工程の出力例を示す。図１９は、特定の母集団Ｕの年度違いの場合であって、各年度Ｑとも８月の推定泌乳量ＥＹが表された場合である。８月の推定泌乳量ＥＹを算出する場合の要因パラメータＰは、８月の特定日の要因パラメータＰであってもよいし、連続した３日若しくは１週間の要因パラメータＰを平均してもよい。

図１９を参照すると、２０１８年から２０２０年にかけ、推定泌乳量ＥＹは全体として高くなっている。特に分娩後日数ｎ_ｍ付近の向上が高くなっている。したがって、各年毎に行なわれた泌乳量改善の効果が出ていることがわかる。

このように推移泌乳量比較ＴＹ工程では、指定された表示期間の推定泌乳量ＥＹを比較することができる。図１９では、異なる年度の同一月の場合の表示例を示したが、同一年度の月違いであってもよい。

また、図２０には、乳牛個体の実乳量ＲＹを同時にプロットした場合の結果について示す。図２０（ａ）は、例えば２０１９年の５月の推定泌乳量ＥＹ（５月）を曲線にした泌乳曲線Ｍ（５月）が表示され、図２０（ｂ）は、２０１９年の８月の推定泌乳量ＥＹ（８月）を曲線にした泌乳曲線Ｍ（８月）が表示されている。なお、図２０（ａ）、（ｂ）共に、他方の月の推定泌乳量を点線で示した。そして、特定の乳牛Ｉ_ｎの実乳量Ｉ_ｎＲＹが表示されている。なお、Ｎ_５０は分娩後日数が５０日の場合を示し、Ｎ_１５０は分娩後日数が１５０日の場合を示している。５月から８月への時間の経過に伴って、個体牛Ｉ_ｎの泌乳量も変化する。

図２０では、この個体牛Ｉ_ｎは５月時点では、推定泌乳量^Ｎ５０ＥＹ（５月）よりも高い泌乳量であったが、８月になると推定泌乳量^Ｎ１５０ＥＹ（８月）並み、若しくはそれより低い値となっている。このように、推移泌乳量比較ＴＹを算出する工程では、特定の乳牛の泌乳量の変化の経緯と、その変化と平均レベルとの相違も知ることができる。

再度図２を参照する。以上のように本発明に係る泌乳量算出システム１は、推定泌乳量ＥＹ、個別推定泌乳量ＩＹ、群推定泌乳量ＧＹ、推定泌乳量比較ＣＹおよび推移泌乳量比較ＴＹを算出することができるが、それぞれの工程の結果を持って、他の工程を引き続き実行できるようにしてもよい。

具体的には、あるメニューの結果を保存しておき、他のメニューを実行している際に、同時に保存しておいた結果を示す若しくは利用することで、メニューの組み合わせは可能となる。

日本国内の牧場で得たデータを用いて本発明に係る泌乳量算出システムを実施した。この牧場は、ある地方の複数の牧場で乳牛数は、およそ５０００頭ほどのデータである。母集団は、２０１９年の月毎のデータである。要因パラメータＰは、ＴＨＩと日照時間である。ＴＨＩは月の平均温度と平均湿度の値を用いた。日照時間は月合計値である。結果を図２１に示す。

図２１（ａ）は、２０１９年の２月、５月、８月、１１月の５か月のＴＨＩと日照時間の値である。対象の牧場の平均実乳量（ｋｇ／日）と分娩後日数ＮをＷＯＯＤ曲線で近似し、そこから得た値を用いて推定泌乳量ＥＹを要因パラメータＰとしてＴＨＩだけを用いた単回帰で求めた結果を図２１（ｂ）に示す。また、要因パラメータＰとして、ＴＨＩに日照時間（月合計）の２つを用いた重回帰によって推定泌乳量ＥＹを求めた結果を図２１（ｃ）に示す。なお、いずれも分娩後日数が５０日の場合である。

要因パラメータＰとしてＴＨＩだけを用いると、標準誤差が０．６５とかなりばらついた。しかし、要因パラメータＰとしてＴＨＩと日照時間を用いると、標準誤差が０．１７となり、かなり実乳量を推定できていることが分かった。

泌乳量がＴＨＩと日照時間で推定できる点は、これまで知られることがなく、今後の泌乳量の増加に対する対応に大きなインパクトを与えるものである。

本発明は、酪農の情報化のために、好適に利用することができる。

１ 泌乳量算出システム
１０ 端末
１２ 本体
１４ メモリ
１６ 牧場
１８ 外部情報

Claims (1)

1. 分娩後日数、実乳量および気候データを含む要因パラメータを牧場データとして蓄積するメモリと、
   前記牧場データの中の一定期間の前記分娩後日数と前記実乳量を補間する補間式を作成する補間式作成部と、
   前記補間式から算出される特定の前記分娩後日数の補間乳量を目的変数とし、前記牧場データから選ばれた少なくとも１つの前記要因パラメータを説明変数として推定泌乳量を算出する回帰式を作成する回帰式作成部と、
   入力された前記要因パラメータの値に対して前記回帰式から推定泌乳量を算出する乳量算出部を含む制御装置を有することを特徴とする泌乳量算出システム。

Images