(57)【要約】【課題】人の骨格モデルの自由度を非常に小さく近似し、荷物の持ち上げ動作に特化した筋張力推定手法は計算コストや直感的な理解のしやすさ等のメリットから労働現場への導入が進められている。多関節拮抗筋を含む筋骨格モデルにおいても一般目的関数を利用した筋張力推定では、筋の共収縮モデルを持たないため関節圧迫力を過小評価する可能性がある。【解決手段】特殊目的関数に係る最適化問題を解いてワイヤ張力を計算し、そのワイヤ張力による関節間力を利用して関節負荷推定を行った。【選択図】図２
2019005344.tif 000002
【特許請求の範囲】
【請求項１】
ワイヤと剛体リンクからなる筋骨格モデルを利用して関節の負荷を推定する方法であって、
前記筋骨格モデルと前記関節がその負荷を受ける運動の情報から前記運動を実現するための前記関節に係る一般化力と、前記運動により生じる前記関節に係る関節間力を逆動力学解析により推定し、
前記一般化力を生じさせるワイヤ張力を求める一般目的関数に、前記関節に係る関節間力の評価項を追加した特殊目的関数に係る最適化問題を解いてワイヤ張力を計算し、
そのワイヤ張力による関節間力を取得することを特徴とする方法。
【請求項２】
請求項１において、前記逆動力学解析により推定された当該関節間力と前記取得されたワイヤ張力による関節間力の合力を計算し、
前記合力を前記関節に係る負荷とすることを特徴とする方法。
【請求項３】
前記関節がその負荷を受ける運動の情報は一般化座標、または一般化座標と速度、加速度、接触力のいずれかまたは全ての情報であって、実際に運動を計測し、または、前記運動のシミュレーションにより得た情報であることを特徴とする請求項１又は請求項２のいずれかに記載の方法。
【請求項４】
前記一般目的関数は、所定の等式拘束条件、不等式拘束条件下で次の数式によることを特徴とする請求項３に記載の方法。ただし、τjは一般化力、その一般化力を生じさせるワイヤ張力をｆ、Wf、Wτは重み行列、JjTはワイヤ張力を一般化力へ変換するための一般化座標qに関するヤコビ行列とする。
【数１１】
2019005344.tif 000013
【請求項５】
前記最適化問題は非線形計画問題であって前記一般目的関数は、所定の等式拘束条件、不等式拘束条件下で次の数式によることを特徴とする請求項３に記載の方法。
【数１２】
2019005344.tif 000014
ここで、次数nとmは3乗または4乗とする。
【請求項６】
前記最適化問題は非線形計画問題であって前記一般目的関数は、所定の等式拘束条件、不等式拘束条件下で次の数式によることを特徴とする請求項３に記載の方法。ここで、次数nは2乗、3乗または4乗とする。
【数１３】
2019005344.tif 000015
【請求項７】
前記関節間力の評価項を加えた特殊目的関数は次の数式によることを特徴とする請求項４乃至６のいずれか１項に記載の方法。
【数１４】
2019005344.tif 000016
【請求項８】
前記合力は次の数式によることを特徴とする請求項７に記載の方法。
【数１５】
2019005344.tif 000017
【請求項９】
請求項１乃至請求項８のいずれか１項において、前記特殊目的関数の各項に係る重みを所望の制約条件下で定めることにより、前記関節に係る負荷の上限を求めることを特徴とする方法。
【請求項１０】
請求項１乃至請求項９のいずれか１項に記載の方法を実行することを特徴とするプログラムおよびプログラムを記録した記憶媒体。
【請求項１１】
動体計測装置と演算装置と表示装置を備えワイヤと剛体リンクからなる筋骨格モデルを利用して関節の負荷を推定するシステムであって、
前記動体計測装置において、実際に前記関節がその負荷を受ける運動を計測して、前記関節に係る一般化座標、または一般化座標と速度、加速度、接触力のいずれかまたは全ての情報を取得し、
前記演算装置において、
前記筋骨格モデルと前記関節がその負荷を受ける運動の前記情報から前記運動を実現するための前記関節に係る一般化力と、前記運動により生じる前記関節に係る関節間力を逆動力学解析により推定し、
前記一般化力を生じさせるワイヤ張力を求める一般目的関数に、前記関節に係る関節間力の評価項を追加した特殊目的関数に係る最適化問題を解いてワイヤ張力を計算し、そのワイヤ張力による関節間力を取得し、
前記逆動力学解析により推定された当該関節間力と前記取得されたワイヤ張力による関節間力の合力を計算し、
前記表示装置において、前記合力を前記関節に係る負荷として表示することを特徴とするシステム。
【請求項１２】
請求項１１に記載のシステムにおいて、前記動体計測装置に代わって、前記演算装置で前記関節に係る一般化座標、または一般化座標と速度、加速度、接触力のいずれかまたは全ての情報を前記運動のシミュレーション解析により求めたことを特徴とするシステム。
【請求項１３】
請求項１１乃至請求項１２のいずれかにおいて、前記特殊目的関数の各項に係る前記重みを所望の制約条件下で定めることにより、前記関節に係る負荷の上限を求めることを特徴とするシステム。
【発明の詳細な説明】
【技術分野】
【０００１】
本発明は筋骨格モデルによる動作解析を用いた関節間力の評価手法に関する。
【背景技術】
【０００２】
人体の内力は外部からの計測が困難であるため、シミュレーションによる負荷推定が研究されている。
【０００３】
労働安全衛生分野では、簡略化した人体モデルを用いて椎間板内圧等を指標として作業姿勢を評価する取り組みが行われてきた（非特許文献１）。
また装着型のデバイスを用いたリアルタイム負荷警告システム等の開発も行われている（非特許文献２）。
【０００４】
これら作業負荷推定手法は、２次元的な動き（腰の屈曲伸展方向）のみ扱い、人の骨格モデルの自由度を非常に小さく近似し、荷物の持ち上げ動作に特化した開発を進めていることが多い。
そのような推定手法は計算コストや直感的な理解のしやすさ等のメリットから労働現場への導入が進められているが、一方で筋の共収縮の効果を考慮していないため関節圧迫力を過小評価する可能性がある。
【０００５】
高橋らは椎間板圧力の実測により背筋のみをモデル化した場合の理論値を計測値が上回ることを確認しており、それが筋活動の影響であると結論付けている（非特許文献３）。
したがって、リスク推定のためにはこの過小評価の問題を解決するような推定方法の開発が必要であると考えられる。
【先行技術文献】
【非特許文献】
【０００６】
【非特許文献１】瀬尾  明彦：“作業関連性筋骨格系障害に関わる計測と評価の手法”  労働科学、81(1)、pp.6-15、2005.
【非特許文献２】Y. Tsuchiya, Y. Matsuo and T. Tanaka: “Estimation of lumbar loadBy 2D reconstruction of spine line using wearable sensor system,” 2014IEEE International Conference on Systems, Man and Cybernetics (SMC), pp.3669-3674, 2014.
【非特許文献３】I. Takahashi, S. Kikuchi, K. Sato and N. Sato: “Mechanical load of the lumbar spine during forward bending motion of the trunk-a biomechanical study,” Spine, 31(1), pp.18-23, 2006.
【非特許文献４】Y. Nakamura, K. Yamane, Y. Fujita, and I. Suzuki: “Somatosensory Computation for Man-Machine Interface from Motion-Capture Data and Musculoskeletal Human Model,” IEEE Trans. on Robotics, 21(1), pp.58-66, 2005.
【発明の概要】
【発明が解決しようとする課題】
【０００７】
筋骨格モデルはヒトの身体を関節によって連結された剛体リンクとしてモデル化し、関節駆動力を発揮するための筋・腱・靭帯などをワイヤとして表現する。
筋骨格モデルを用いたシミュレーションでは計測した動作や外力を入力とし、逆動力学解析によりその運動を実現するような関節トルクやワイヤの張力を求める。
【０００８】
ただし、ヒトの筋は自由度に対して冗長であるため、ワイヤ張力を一意に決定することは出来ない。
そのため筋張力の２乗和あるいは３乗和を目的関数とした最適化問題により筋張力を求めるなど筋負担を平均化・最小化することが一般的である。
【０００９】
しかし、この方法では動作に必要な最低限の筋活動が解として得られるため身体負担としては過小評価となる可能性がある。
【課題を解決するための手段】
【００１０】
本発明は、 多関節拮抗筋を含む筋骨格モデルによる関節負荷推定方法を提供できる。
（１）
ワイヤと剛体リンクからなる筋骨格モデルを利用して関節の負荷を推定する方法であって、
前記筋骨格モデルと前記関節がその負荷を受ける運動の情報から前記運動を実現するための前記関節に係る一般化力と、前記運動により生じる前記関節に係る関節間力を逆動力学解析により推定し、
前記一般化力を生じさせるワイヤ張力を求める一般目的関数に、前記関節に係る関節間力の評価項を追加した特殊目的関数に係る最適化問題を解いてワイヤ張力を計算し、
そのワイヤ張力による関節間力を取得することを特徴とする方法。
（２）
（１）において、前記逆動力学解析により推定された当該関節間力と前記取得されたワイヤ張力による関節間力の合力を計算し、
前記合力を前記関節に係る負荷とすることを特徴とする方法。
（３）
前記関節がその負荷を受ける運動の情報は一般化座標、または一般化座標と速度、加速度、接触力のいずれかまたは全ての情報であって、実際に運動を計測し、または、前記運動のシミュレーションにより得た情報であることを特徴とする（１）又は（２）のいずれかに記載の方法。
【００１１】
（４）
前記一般目的関数は、所定の等式拘束条件、不等式拘束条件下で数式（５）によることを特徴とする（３）に記載の方法。ただし、τjは一般化力、その一般化力を生じさせるワイヤ張力をｆ、Wf、Wτは重み行列、JjTはワイヤ張力を一般化力へ変換するための一般化座標qに関するヤコビ行列とする。
（５）
前記最適化問題は非線形計画問題であって前記一般目的関数は、所定の等式拘束条件、不等式拘束条件下で数式（９）によることを特徴とする（３）に記載の方法。ここで、次数nとmは3乗または4乗とする。
（６）
前記最適化問題は非線形計画問題であって前記一般目的関数は、所定の等式拘束条件、不等式拘束条件下で数式（１０）によることを特徴とする（３）に記載の方法。ここで、次数nは2乗、3乗または4乗とする。
（７）
前記関節間力の評価項を加えた特殊目的関数は数式（６）によることを特徴とする（４）乃至（６）のいずれかに記載の方法。
（８）
【００１２】
前記合力は数式（７）によることを特徴とする（７）に記載の方法。
（９）
（１）乃至（８）のいずれか１項において、前記特殊目的関数の各項に係る重みを所望の制約条件下で定めることにより、前記関節に係る負荷の上限を求めることを特徴とする方法。
（１０）
（１）乃至（９）のいずれか１項に記載の方法を実行することを特徴とするプログラムおよびプログラムを記録した記憶媒体。
【００１３】
（１１）
動体計測装置と演算装置と表示装置を備えワイヤと剛体リンクからなる筋骨格モデルを利用して関節の負荷を推定するシステムであって、
前記動体計測装置において、実際に前記関節がその負荷を受ける運動を計測して、前記関節に係る一般化座標、または一般化座標と速度、加速度、接触力のいずれかまたは全ての情報を取得し、
前記演算装置において、
前記筋骨格モデルと前記関節がその負荷を受ける運動の前記情報から前記運動を実現するための前記関節に係る一般化力と、前記運動により生じる前記関節に係る関節間力を逆動力学解析により推定し、
前記一般化力を生じさせるワイヤ張力を求める一般目的関数に、前記関節に係る関節間力の評価項を追加した特殊目的関数に係る最適化問題を解いてワイヤ張力を計算し、そのワイヤ張力による関節間力を取得し、
前記逆動力学解析により推定された当該関節間力と前記取得されたワイヤ張力による関節間力の合力を計算し、
前記表示装置において、前記合力を前記関節に係る負荷として表示することを特徴とするシステム。
（１２）
（１１）に記載のシステムにおいて、前記動体計測装置に代わって、前記演算装置で前記関節に係る一般化座標、または一般化座標と速度、加速度、接触力のいずれかまたは全ての情報を前記運動のシミュレーション解析により求めたことを特徴とするシステム。
（１３）
（１１）乃至（１２）のいずれかにおいて、前記特殊目的関数の各項に係る前記重みを所望の制約条件下で定めることにより、前記関節に係る負荷の上限を求めることを特徴とするシステム。
【発明の効果】
【００１４】
本発明の利点は、最大の関節間力を解析することにより、計測された作業動作の潜在的なリスクを評価できることである。
重量物持ち上げ作業時の椎間板圧縮力の負荷評価に有効であることがわかった。
【図面の簡単な説明】
【００１５】
【図１】本発明の解析対象とした一連の荷物の持ち上げ動作を表す図である。
【図２】腰椎に影響を与える体幹周りの筋と関節のモデルの図である。左図は背面、右図は前面からみた図である。
【図３】解析結果の図である。上図は体幹の鉛直方向からの傾き、下図は腰椎圧縮力の推定値を表す。
【図４】動作開始から９秒後のフレームの体幹の筋張力算出結果を表す図であって、(a)高リスク推定結果、(b)関節間力に対する条件無し、(c)低リスク推定結果を表す図である。
【図５】本発明のフローチャート図を表す。
【発明を実施するための形態】
【００１６】
本発明は、筋骨格モデルを用いた関節間力の推定手法を提案し、実施例として、持ち上げ動作中の椎間板への負担推定方法を説明する。
最初に関節間力を考慮しない一般的なワイヤ張力推定手法について述べ、 それを関節間力評価に拡張した本発明について述べる。
以下は椎間板圧縮力の解析を例にとり説明するが本発明はそれに限定されるものではない。
【００１７】
図５に本発明のフローチャート図を示す。
まず本発明で扱う筋骨格モデルは、人の骨格は適当な細かさで分割して剛体リンクとし、剛体リンク間を回転・球面関節により結合することにより骨格系をモデル化する。
関節駆動力を発揮するための筋・腱・靭帯などは骨格上に付着点を持つワイヤとして表現する（Ｓ１０）。
骨格上に付着点を持つワイヤとして近似された筋・腱・靭帯は、収縮方向に力を発揮する。
人体リンクモデルの運動方程式は式（１）のように表される。
【００１８】
【数１】
2019005344.tif 000003
【００１９】
ただし、JjTはワイヤ張力を一般化力へ変換するためのqに関するヤコビ行列、JEiTは接触点Eiにおける外力を一般化力へ変換するためのqに関するヤコビ行列とする。
また、qは一般化座標、q（上に１点）はその１次微分、q（上に２点）はその２次微分とし、Mは慣性行列、bはコリオリ力・遠心力・重力項、fはワイヤ張力、fEiは接触点Eiに働く外力とする。
【００２０】
まず、一般化座標、速度、加速度、接触点の位置と接触力は、全身の動作と外力を何らかの計測方法により計測する、または解析者が作成・推定することにより既知とする（Ｓ２０）。
この運動を実現するための一般化力τjを次式のように逆動力学解析により計算する（Ｓ３０）。
本研究では中村らによる筋骨格モデルの動力学計算ライブラリを用いた（非特許文献４）。
【００２１】
【数２】
2019005344.tif 000004
【００２２】
さらに、この一般化力τjを生じさせるワイヤ張力、f=[f1 f2 ... fNw]を求める。
【００２３】
【数３】
2019005344.tif 000005
【００２４】
ただしワイヤは収縮方向にのみ力を発揮するため各要素の値はゼロ以下という制約条件を持つ。
特にfj が筋による張力を表すときは次式（４）が条件となる。
【００２５】
【数４】
2019005344.tif 000006
【００２６】
ここで、fmax は各筋の最大筋張力であり解剖学的知見に基づき計算される値である。
【００２７】
ワイヤの数は一般化力に対して冗長であるため、fは一意に定まらない。そのため生体力学的に意味を持つ指標を導入し最適化計算を行う。
本手法では式（５）のような一般目的関数を定め、式（４）を不等式拘束条件とした二次計画問題を解くことによりワイヤ張力を求める。
【００２８】
【数５】
2019005344.tif 000007
【００２９】
ここで、Wf、Wτは重み行列であり、正の値を対角成分に持つ対角行列である。
第１項はワイヤ張力の二乗和を最小化するための項、第２項は式（３）の誤差の評価項である。
第２項は理想的には等式制約条件とすることができるが、実際には運動や外力の計測誤差やモデル化誤差を含んでいるため評価関数に組み込むことにより条件を緩和させた。
【００３０】
以上、関節間力を考慮しない一般的なワイヤ張力推定手法について述べた。
【００３１】
次に、それを関節間力評価に拡張した本発明について述べる。
本発明では、式（５）の目的関数に代わって関節間力評価を行う項を追加した次式（６）の特殊目的関数により解析を行う（Ｓ４０）。
【００３２】
【数６】
2019005344.tif 000008
【００３３】
一般的なワイヤ張力推定手法を表す式（５）では、関節間に生じる力のうち、運動に関与する一般化力のみを用いている。
【００３４】
しかし、実際には与えられた運動（q、q（上に１点）、q（上に２点））と外力fEiから、各関節に3軸方向の並進力Fと、各軸周りのモーメントNが逆動力学解析により計算される。
【００３５】
この力とモーメントを受ける軸が駆動軸の場合には関節に運動を生じさせるが、固定軸の場合にはつり合う拘束力が生じるため運動としては現れない。
最終的な関節間力は自重や外力による関節間力τcとワイヤ張力の影響による関節間力の合力、
【００３６】
【数７】
2019005344.tif 000009
【００３７】
である（Ｓ５０）。
【００３８】
この関節間力のうち、椎間板に働く脊柱方向の力が椎間板圧縮力である（Ｓ６０）。
ステップＳ４０の式（６）における目的関数の重みWcが正である場合は対応する関節間力を減少させ、負にした場合は関節間力を増加させる。
【００３９】
つまり、椎間板圧縮力に対応する重みを負にした場合、最も椎間板圧縮力が大きくなるような筋張力の解とそのときの圧縮力を求めることが出来る。
この重みが負のときを関節への負荷が大きくなる高リスク解析、正のときを低リスク解析と定義すればよい。
すなわち、式（６）の２次計画問題の重みの設計により、高リスク・低リスクなど必要な解析を行う。
必ずしも高リスク解析・低リスク解析を常にどちらも行う必要はない。
【００４０】
この重みを大きくするほど一般化力の誤差が増大するため、解析ごとに許容誤差を定めてステップＳ４０の式（６）における重み行列を適切に設定する必要がある。
【実施例１】
【００４１】
提案手法を持ち上げ動作中の椎間板負荷の解析に適用した例を示す。
図1に解析対象とした荷物の持ち上げ動作を示す。
動作は光学式モーションキャプチャシステムを用いて計測した。
【００４２】
使用した人体モデルは体重54.4kg、身長1.68m、47自由度、ワイヤ総数は360本である。両手にそれぞれ2.5kgの質量を加え、合わせて5kgの手先負荷を与えた。
【００４３】
図２に腰椎に影響を与える体幹周りの筋のモデルを示す。
体幹のリンクは腰仙椎境界（第５腰椎L5と第1仙骨S1間）と、胸腰椎境界（第12胸椎T12と第1腰椎L1の間）で分割し、骨盤、腰椎、胸椎の３つのリンクとした。
【００４４】
目的関数(6)の重みは実験的に次のように設定した。
Wτの対角要素はwτ=1。
Wfの対角要素は、
【００４５】
【数８】
2019005344.tif 000010
【００４６】
ただし、fjmax は生理学的断面積などの知見を元に定める各筋の最大張力である。
ここではwf=1とした。
【００４７】
Wcの対角要素wcは、L5/S1間の椎間板の圧縮方向のみに値を与え、他の項は0とした。
高リスク解析（圧縮力最大化）のときwc=-3×10-5、低リスク解析（圧縮力最小化）のときwc=3×10-5とした。
【００４８】
比較のため、関節間力を考慮しない式（５）の目的関数による筋張力推定も行い、これを通常の筋張力解析結果とする。
各項目の重みはwτ=1,wf =20とした。
【００４９】
図３に椎間板圧縮力の解析結果を示す。
【００５０】
図３の上段のグラフは、体幹の前傾角の時系列を表す。
図３の下段のグラフは、各時系列におけるL5/S1間の椎間板圧縮力の推移を表す曲線であり、上から順に高リスク推定結果、通常の筋張力解析の結果、低リスク推定結果を表す曲線を表している。
【００５１】
通常の筋張力推定によって求めた椎間板圧縮力に対し、動作全体に渡って低リスク条件では値が小さく、高リスク条件では値が大きくなっていることが確認できる。
【００５２】
図４に動作開始から９秒後のフレームの体幹の筋張力算出結果を示す。
【００５３】
ここで、仮に脊柱をひとつの関節でモデル化し1対の拮抗筋によって駆動すると仮定した場合、椎間板圧縮力を最大化しようとすると主動筋または拮抗筋のどちらかが必ず最大張力となる。
【００５４】
つまりそのようなモデル化を行った場合、筋が自然長に近く力を発揮しやすい直立姿勢で最も負荷が大きくなる。
【００５５】
しかし、今回のように複数の関節をまたぐ筋によるモデル化を行う場合、その筋張力は複数の関節の力のつりあいを同時に維持する必要があるため、必ずしも全ての姿勢で最大張力にはならない。
これによって、より自然な共収縮の状況を再現できていると考えられる。
【００５６】
図４を見ると、図４(b)の関節間力を考慮しない条件では背筋に最大筋張力fjmaxに応じた筋張力が分配されているのに対し、図４(c)の低リスク条件では椎間板の圧縮方向への影響の大きい筋の出力が減少している。
【００５７】
また、図４(a)の高リスク条件では椎間板の圧縮方向に沿った筋の活性度が上がると共に、腹筋の筋活動も上昇している。
実際の椎間板への圧縮力はこの低リスク推定結果と高リスク推定結果の間に存在すると考えられる。
【００５８】
以上まとめると、提案手法により計算された椎間板圧縮力は一般的手法によって計算された負荷を上回っており、積極的な共収縮を再現した解析結果が得られることが確認された。
【実施例２】
【００５９】
実施例１では、二次計画問題によりワイヤ張力を求めたが、一般目的関数を次式のように設定し、非線形計画問題によりワイヤ張力を推定してもよい。
【００６０】
【数９】
2019005344.tif 000011
【００６１】
ここで、次数nとmは3乗または4乗である。
また、式（９）の第２項を等式制約条件とし、次の一般目的関数によりワイヤ張力を推定してもよい。
【００６２】
【数１０】
2019005344.tif 000012
【００６３】
ここで、次数nとmは3乗または4乗である。
一般目的関数は式（５）、式（９）、式（１０）以外の、ワイヤ張力を導き出せるどのような目的関数であってもよい。
【００６４】
関節間力を取り入れた特殊目的関数はいずれの一般目的関数を用いても数式（６）のような関数としてよい。
【実施例３】
【００６５】
実施例１では、動作は光学式モーションキャプチャシステムを用いて計測したが、リアルタイム動作計測技術と本手法を組み合わせれば、実際の作業中にこのような指標に基づきリスクが高いと判断される姿勢を取ったときに警告を発することが可能となる。
【００６６】
特殊目的関数の重みの比により一般化力の誤差が変化するため、その閾値を定め、誤差が閾値を超えない範囲で対象の関節間力に係る重みを最小とするような重み付けをすることにより、対象とする関節に想定される最大負荷を求めることが可能となる。
【産業上の利用可能性】
【００６７】
本発明は特に運動時の関節へのリスク推定を目的とした解析において有用であると考えられる。
例えば、衝撃を予測して身体を強張らせているときとリラックスしているときでは、見かけ上は同じ姿勢・運動をしているとしても内力は異なっている。
また、プロスポーツ選手とアマチュア選手や、作業の熟練者や初心者など個人差による力の入れ方の違い等も存在する。
それらの解析を本手法により扱うことが出来る。
【符号の説明】
【００６８】
１  胸最長筋
２  腰腸肋筋
３  腰方形筋
４  腹直筋
５  外腹斜筋
６  内腹斜筋
７  T12/L1関節
８  L5/S1関節
【図１】
2019005344.tif 000018
【図２】
2019005344.tif 000019
【図３】
2019005344.tif 000020
【図４】
2019005344.tif 000021
【図５】
2019005344.tif 000022