(57)【要約】【課題】不良発生メカニズムにかかわる統計モデルを同定すること。【解決手段】各不良モードの不良発生メカニズムに関する仮説を複数個想定した上で、仮説のそれぞれについて、複数の状態変数から複数の出力変数へのデータ分布の関数型を、各不良モードに対応した所定の故障物理モデルを考慮して仮定し、複数の状態変数に関する確率変数データ集合と離散化した結合確率分布および確率密度が０でない確率変数データ集合とを、あらかじめ与えられたまたは数値実験点を拡充中のＣＡＥ解析結果データに基づき応答曲面モデルおよびモンテカルロシミュレーションを用いて算出し、複数の入力変数および前記複数の出力変数についてのあらかじめ得られたまたはモニタリング中の観測データ集合と、第２のステップにおいて算出された複数の状態変数の結合確率分布とから、仮定されたデータ分布の関数形に応じたパラメータを推定する。【選択図】図３
2007108843.tif 000002
【特許請求の範囲】
【請求項１】
不良発生にかかわる複数の設計変数（不良要因）をそれぞれ入力変数とし複数の不良モードにかかわる測定可能量をそれぞれ出力変数とし、かつ、前記複数の入力変数と前記複数の出力変数との間に複数の潜在的な状態変数が介在するとした階層構造の統計モデルにおいて、各不良モードの不良発生メカニズムに関する仮説を複数個想定した上で、前記仮説のそれぞれについて、前記複数の状態変数から前記複数の出力変数へのデータ分布の関数型を、各不良モードに対応した所定の故障物理モデルを考慮して仮定する第１のステップと、
前記複数の状態変数に関する確率変数データ集合と離散化した結合確率分布および確率密度が０でない確率変数データ集合とを、あらかじめ与えられたまたは数値実験点を拡充中のＣＡＥ解析結果データに基づき応答曲面モデルおよびモンテカルロシミュレーションを用いて算出する第２のステップと、
前記複数の入力変数および前記複数の出力変数についてのあらかじめ得られたまたはモニタリング中の観測データ集合と、前記第２のステップにおいて算出された複数の状態変数の結合確率分布とから、前記仮定されたデータ分布の関数形に応じたパラメータを推定する第３のステップと
を具備することを特徴とする半導体装置設計支援方法。
【請求項２】
前記第３のステップによって得られた階層構造の統計モデルから、適合度の高い統計モデルを選定し、該選定された統計モデルを用いて、前記不良発生メカニズムの仮設に関する不良要因（入力変数）を同定する第４のステップをさらに具備することを特徴とする請求項１記載の半導体装置設計支援方法。
【請求項３】
前記第１のステップが、前記データ分布として比例ハザードモデルを考慮した階層ベイズモデルを設定し、
前記第３のステップが、ＭＡＰ推定法、最尤法、ＥＭ法、ブートストラップ法を利用してなされることを特徴とする請求項１記載の半導体装置設計支援方法。
【請求項４】
不良発生にかかわる複数の設計変数（不良要因）をそれぞれ入力変数とし複数の不良モードにかかわる測定可能量をそれぞれ出力変数とし、かつ、前記複数の入力変数と前記複数の出力変数との間に複数の潜在的な状態変数が介在するとした階層構造の統計モデルにおいて、各不良モードの不良発生メカニズムに関する仮説を複数個想定した上で、前記仮説のそれぞれについて、前記複数の状態変数から前記複数の出力変数へのデータ分布の関数型を、各不良モードに対応した所定の故障物理モデルを考慮して仮定する第１の手段と、
前記複数の状態変数に関する確率変数データ集合と離散化した結合確率分布および確率密度が０でない確率変数データ集合とを、あらかじめ与えられたまたは数値実験点を拡充中のＣＡＥ解析結果データに基づき応答曲面モデルおよびモンテカルロシミュレーションを用いて算出する第２の手段と、
前記複数の入力変数および前記複数の出力変数についてのあらかじめ得られたまたはモニタリング中の観測データ集合と、前記第２のステップにおいて算出された複数の状態変数の結合確率分布とから、前記仮定されたデータ分布の関数形に応じたパラメータを推定する第３の手段と
を具備することを特徴とする半導体装置設計支援システム。
【請求項５】
不良発生にかかわる複数の設計変数（不良要因）をそれぞれ入力変数とし複数の不良モードにかかわる測定可能量をそれぞれ出力変数とし、かつ、前記複数の入力変数と前記複数の出力変数との間に複数の潜在的な状態変数が介在するとした階層構造の統計モデルにおいて、各不良モードの不良発生メカニズムに関する仮説を複数個想定した上で、前記仮説のそれぞれについて、前記複数の状態変数から前記複数の出力変数へのデータ分布の関数型を、各不良モードに対応した所定の故障物理モデルを考慮して仮定する第１のステップと、
前記複数の状態変数に関する確率変数データ集合と離散化した結合確率分布および確率密度が０でない確率変数データ集合とを、あらかじめ与えられたまたは数値実験点を拡充中のＣＡＥ解析結果データに基づき応答曲面モデルおよびモンテカルロシミュレーションを用いて算出する第２のステップと、
前記複数の入力変数および前記複数の出力変数についてのあらかじめ得られたまたはモニタリング中の観測データ集合と、前記第２のステップにおいて算出された複数の状態変数の結合確率分布とから、前記仮定されたデータ分布の関数形に応じたパラメータを推定する第３のステップと
をコンピュータに実行させるための半導体装置設計支援プログラム。
【発明の詳細な説明】
【技術分野】
【０００１】
本発明は、半導体装置設計支援方法、半導体装置設計支援システム、および半導体装置設計支援プログラムに係り、特に、不良発生メカニズムに関する仮説検証を行うのに有用な半導体装置設計支援方法、半導体装置設計支援システム、および半導体装置設計支援プログラムに関する。
【背景技術】
【０００２】
設計段階において信頼性に関するデザインレビューを強化するためには、現象を再現できる物理モデルに基づくＣＡＥ（Computer Aided Engineering）や、材料、製造、検査・試験、品質に関するデータベースを活用し、設計に内在する不確定性を統計・確率論的手法により考慮した上で、信頼性設計を行う必要がある。
【０００３】
半導体装置においては、要素解析技術や検査・評価技術の進歩および、材料・製造プロセスにおける不良モードに関する知識の蓄積などによって、シミュレーション技術を利用した開発が行なわれている。主に、製造プロセスや試験・検査および市場条件における複雑な現象を再現することにより不良メカニズムを明らかにするため、あるいは、シミュレーション結果と不良発生基準を照らし合わせて材料／構造／プロセスの最適化を行うためなどの目的で活用されている。
【０００４】
半導体装置開発においては、多くの要因が複雑に絡み合いながら複数の不良モードが存在する場合も多く、すべての組み合わせを実験により忠実に再現することは容易ではない。また、例えば半導体パッケージの素子内部や微細接合部のような微小領域では、温度・応力・ひずみなどの物理量を実際に測定して評価することが困難な場合が多い。そのため、シミュレーション技術に基づく問題解決への期待は大きい。最近では、設計に内在するばらつきに対してロバストな設計を行うため、あるいは安全係数（設計マージン）を合理的に設定するためなどの信頼性設計にもシミュレーション技術が適用されつつある。
【０００５】
しかし、ＣＡＥモデルには、特に材料の異方性特性や非弾性特性などのモデル化における構成方程式や、対象とする数値計算モデルの時・空間的な境界条件において不確定性が内在する。バルク材やＴＥＧ（Test Element Group）による個別的な実験データや、専門知識や経験による判断のみでは、限られた費用と期間の中で、各不良モードに対する現象発生のメカニズムを明確化できずに、解析モデルの不確定性を同定しきれない場合も少なくない。
【０００６】
また、調達・製造のグローバル化の中で、不良モードは多様化し、材料特性や製造プロセスのばらつきは大きくなる傾向にある。さらに、半導体装置の微細化・信号伝送高速化・大容量化のための技術の革新や新たなコンセプトに基づく半導体装置の登場は、時間的に不連続な形で起こる。例えば、新しい材料や構造の開発あるいは新たな環境下や異なる負荷想定の下での使用などである。このとき、必要十分な信頼性スペックを、過去の経験を外挿して決定することが困難となり大きな不確定性が残る。
【０００７】
なお、本願に関連する技術を記載した文献には下記のものがある。この文献の開示では、複数の設計変数と複数の出力変数との間に複数の潜在的な状態変数が介在するとした階層構造の考え方はなく、当然その統計（確率）モデルやこれを同定することについても言及はない。
【特許文献１】特開平６−１９５４０２号公報
【発明の開示】
【発明が解決しようとする課題】
【０００８】
上記で述べた状況から、半導体装置設計では、物理現象モデルに基づく不良現象解析・評価技術の深耕と蓄積に加えて、多くの不良メカニズムに関する仮説を、不良要因に関する不確定性が内在する中で合理的に検証する仕組みを構築することが求められる。また、材料や製造および品質に関する全体的なデータから、故障物理モデルに立脚した統計モデルやＣＡＥを活用することにより、不良メカニズムの同定を支援するための新たな手法が求められる。
【０００９】
本発明は、半導体装置設計支援方法、半導体装置設計支援システム、および半導体装置設計支援プログラムにおいて、ＣＡＥデータと実際の測定で蓄積されるデータとをトータルなデータとして用いて、不良発生メカニズムにかかわる統計モデルを同定することが可能な半導体装置設計支援方法、半導体装置設計支援システム、および半導体装置設計支援プログラムを提供することを目的とする。
【課題を解決するための手段】
【００１０】
本発明の一態様に係る半導体装置設計支援方法は、不良発生にかかわる複数の設計変数（不良要因）をそれぞれ入力変数とし複数の不良モードにかかわる測定可能量をそれぞれ出力変数とし、かつ、前記複数の入力変数と前記複数の出力変数との間に複数の潜在的な状態変数が介在するとした階層構造の統計モデルにおいて、各不良モードの不良発生メカニズムに関する仮説を複数個想定した上で、前記仮説のそれぞれについて、前記複数の状態変数から前記複数の出力変数へのデータ分布の関数型を、各不良モードに対応した所定の故障物理モデルを考慮して仮定する第１のステップと、前記複数の状態変数に関する確率変数データ集合と離散化した結合確率分布および確率密度が０でない確率変数データ集合とを、あらかじめ与えられたまたは数値実験点を拡充中のＣＡＥ解析結果データに基づき応答曲面モデルおよびモンテカルロシミュレーションを用いて算出する第２のステップと、前記複数の入力変数および前記複数の出力変数についてのあらかじめ得られたまたはモニタリング中の観測データ集合と、前記第２のステップにおいて算出された複数の状態変数の結合確率分布とから、前記仮定されたデータ分布の関数形に応じたパラメータを推定する第３のステップとを具備する。
【００１１】
また、本発明の一態様に係る半導体装置設計支援装置は、不良発生にかかわる複数の設計変数（不良要因）をそれぞれ入力変数とし複数の不良モードにかかわる測定可能量をそれぞれ出力変数とし、かつ、前記複数の入力変数と前記複数の出力変数との間に複数の潜在的な状態変数が介在するとした階層構造の統計モデルにおいて、各不良モードの不良発生メカニズムに関する仮説を複数個想定した上で、前記仮説のそれぞれについて、前記複数の状態変数から前記複数の出力変数へのデータ分布の関数型を、各不良モードに対応した所定の故障物理モデルを考慮して仮定する第１の手段と、前記複数の状態変数に関する確率変数データ集合と離散化した結合確率分布および確率密度が０でない確率変数データ集合とを、あらかじめ与えられたまたは数値実験点を拡充中のＣＡＥ解析結果データに基づき応答曲面モデルおよびモンテカルロシミュレーションを用いて算出する第２の手段と、前記複数の入力変数および前記複数の出力変数についてのあらかじめ得られたまたはモニタリング中の観測データ集合と、前記第２のステップにおいて算出された複数の状態変数の結合確率分布とから、前記仮定されたデータ分布の関数形に応じたパラメータを推定する第３の手段とを具備する。
【００１２】
また、本発明の一態様に係る半導体装置設計支援プログラムは、不良発生にかかわる複数の設計変数（不良要因）をそれぞれ入力変数とし複数の不良モードにかかわる測定可能量をそれぞれ出力変数とし、かつ、前記複数の入力変数と前記複数の出力変数との間に複数の潜在的な状態変数が介在するとした階層構造の統計モデルにおいて、各不良モードの不良発生メカニズムに関する仮説を複数個想定した上で、前記仮説のそれぞれについて、前記複数の状態変数から前記複数の出力変数へのデータ分布の関数型を、各不良モードに対応した所定の故障物理モデルを考慮して仮定する第１のステップと、前記複数の状態変数に関する確率変数データ集合と離散化した結合確率分布および確率密度が０でない確率変数データ集合とを、あらかじめ与えられたまたは数値実験点を拡充中のＣＡＥ解析結果データに基づき応答曲面モデルおよびモンテカルロシミュレーションを用いて算出する第２のステップと、前記複数の入力変数および前記複数の出力変数についてのあらかじめ得られたまたはモニタリング中の観測データ集合と、前記第２のステップにおいて算出された複数の状態変数の結合確率分布とから、前記仮定されたデータ分布の関数形に応じたパラメータを推定する第３のステップとをコンピュータに実行させる。
【発明の効果】
【００１３】
本発明によれば、半導体装置設計支援方法、半導体装置設計支援システム、および半導体装置設計支援プログラムにおいて、ＣＡＥデータと実際の測定で蓄積されるデータとをトータルなデータとして用いて、不良発生メカニズムにかかわる統計モデルを同定することができる。
【発明を実施するための最良の形態】
【００１４】
以下では本発明の実施形態を図面を参照しながら説明する。まず、図１は、複数の設計変数を入力変数とし複数の不良モードにかかわる測定可能量を出力変数とし、かつ、それらの間に複数の潜在的な状態変数が介在するとした階層構造の統計モデルを説明するための概観図である。以下で説明する実施形態は、この図示の階層構造の想定が前提となる。
【００１５】
図１に示すように、設計対象となる半導体装置において、観測され得る量（測定可能量）は各不良モードｙ１…に対応する。ｙ１…は出力変数１３でもある。各不良モードｙ１…の要因ｘ１…は複数あり、図示するように、負荷ばらつきなどのばらつき要因により、それぞればらつきを有している。要因ｘ１…は入力変数（設計変数）１１でもある。出力変数１３および入力変数１１は、観測データ１５により具体的な値が与えられ得る。
【００１６】
入力変数１１と出力変数１３との間には、各不良モードｙ１…を抑制するために知得したいものの実際には測定できない量（状態量ｚ１…）が介在しているものとする。状態量ｚ１…は状態変数１２である。状態変数１２は、設計、生産などのためのシミュレーション技術（いわゆるＣＡＥ）で得られるＣＡＥデータ１４により具体的な値が与えられ得る。状態量ｚ１…から各不良モードｙ１…への関係には、故障物理モデル１６が適用され得る。
【００１７】
本発明の一実施形態に係る半導体装置設計支援方法をまず概略的に説明する。複数の不良モードにおける観測データ１５と、不良要因に関するＣＡＥを活用したパラメータサーベイ結果（ＣＡＥデータ１４）とから、観測データを説明し不良メカニズムを同定するための統計モデルを作成するために、階層ベイズモデルと応答曲面法を活用する。ここで階層ベイズモデルにおける潜在変数（状態量）に関する分布は、ＣＡＥに基づく応答曲面モデルからモンテカルロシミュレーションにより求めた確率分布を用いる。また、データ分布モデルのパラメータ推定には最大事後確率推定（ＭＡＰ（Maximization A Posterior）推定）やＥＭ（Expectation Maximization）法などの統計的推定法を用い、モデルのロバスト性評価についてはリサンプリング法を活用する。これにより、観測データとＣＡＥ結果を融合したかたちで、不良メカニズムを表現する統計モデルの推定が可能となる。
【００１８】
以上のデータ処理の概略フローは以下である。
１）故障物理（図２を参照）を考慮して観測データのデータ分布モデルｐ（ｙ｜ｚ；α）を設定
ここで比例ハザードモデルを用いる場合：
【数１】
2007108843.tif 000003
ただし、
【数２】
2007108843.tif 000004
２）ＣＡＥと応答曲面モデルとモンテカルロシミュレーションを活用し、潜在変数ｚに関する確率分布ｐ（ｚ｜ｗ）を算出
１．応答曲面モデルとモンテカルロシミュレーションにより状態量の確率分布を計算
２．分割範囲と分割数を指定後、確率分布の離散化
３．確率密度が０でないデータ集合ｚを抽出
３）事前分布のハイパーパラメータα、ｗを選定し、データ分布のモデルパラメータθを統計的に推定
ここでパラメータθの選択や識別には情報量基準に基づく方法を活用する。
また本推定法においてブートストラップ法を活用することもできる。
４） ベイズ法により事後分布ｐ（ｚ｜ｙ；α、ｗ）を算出
【数３】
2007108843.tif 000005
５）推定した統計モデルの安定性を評価（推定条件に摂動を与えることにより評価）
６）得られた事後分布から不良発生仮説を同定（本実施形態ではカルマンフィルタと多重仮説検定を活用）
【００１９】
図３は、本発明の一実施形態に係る半導体装置設計支援方法を示すやや詳細なフロー図である。以下、図３に基づき説明する。
【００２０】
対象とする複合あるいは単一の不良現象において、まず、不良発生メカニズム（仮説）を想定することが前提となる（ステップ３０、不良現象データベース２１）。考えられる仮説の想定は、各不良モードｙ１…と不良要因ｘ１…の洗い出しの後に行う。ただし、仮説の想定においては、多くの不良要因ｘ１…に関する組合せや不確定性の中で考えられる仮説だけではなく、ＣＡＥ解析モデルにおける物理現象モデリング、材料モデリング、境界条件、負荷条件などの不確定性の中で考えられる要因の組合せとして得られる仮説も含む。
【００２１】
対象とする各不良モードｙ１…について、各仮説ごとに、対象とする不良要因ｘ１と潜在変数（状態量ｚ１…）を決定する（ステップ５１、５３）。また、各不良モードｙ１…ごとに、対象とする観測データの決定も行う（ステップ５２）。これらの不良モードｙ１や不良要因ｘ１および観測データの決定は、観測データベース４４を参照することにより行う。ここで潜在変数とは、実際には精度の良い実測を容易に実施することが困難であるが、不良現象を左右する重要な状態量（変形量、ひずみ範囲、曲率、応力など）を意味している。
【００２２】
次に、各不良モードｙ１…に対するデータ分布を、対象とする不良モードｙ１…における故障物理モデルを考慮した上で決定する（ステップ５４、データ分布・故障物理モデルベース７１）。このデータ分布の統計モデル（確率モデル）には、複数あるいは単一の、観測変数ｙ、潜在変数（状態量）ｚおよびモデルパラメータθが含まれる。故障物理モデルには、一般に、大別してストレス−ストレングスモデルと確率モデルとがあり、前者の例として破壊機構モデル、速度過程モデル、累積損傷モデルなどがあり、後者の例として故障率モデル、多状態確率過程モデル（マルコフモデル、マルチステートモデル）、最弱リンクモデル（ワイブルモデル）、束モデル（バンドルモデル、ローゼンモデル）などがある（図２を参照）。
【００２３】
データ分布の統計モデルの想定例として、比例ハザードモデルを用いる場合について以下で説明する。ただし、比例ハザードモデルだけではなくマルコフモデルなどの他の故障物理に立脚した統計モデルにも適用できるため、各不良現象の故障物理に即したモデルを選択するのが望ましい。数多くの不良現象に関わる不良要因ｘ１…と不良データの間に明確で再現性のある因果関係が見出せない場合においても、不良現象を支配する状態量（潜在変数）ｚ１…と不良データの関係には、ある種の確率的依存性が見出される場合が多い。
【００２４】
この確率的依存性を、個々の変数が複合的なリスクとして影響を与える競合リスク問題として扱い、各変数の相対的重要度を観測データに基づき判断するための統計的推定モデルとして比例ハザードモデルが存在する。比例ハザードモデルは、不良現象に影響を与える変数の効果を定量的に評価し、その重要度を決定することができる。故障率μは、基線故障率μ0に対して複数の変数の効果を多重的に受けることを表現するために以下のように表す。ただし、変数間に交互作用がある場合には、交互作用項を以下のexp関数の引数の中に含めてもよい。
【数４】
2007108843.tif 000006
このとき、不良発生確率は以下のようにかける。
【数５】
2007108843.tif 000007
ここで、Ｂは統計モデルのパラメータであり、観測データをもとに推定される。また、Ｒ0とＲ（ｔ｜ｚ）は以下の式であり、故障物理を考慮した上で分布型を決める。
【数６】
2007108843.tif 000008
【数７】
2007108843.tif 000009
【００２５】
この比例ハザードモデルには、Ｒ0に含まれる未知パラメータと、Ｂｉの未知パラメータが内在する。また、不良データが寿命でない場合には、tではなく強度データyなどの物理量となる。
【００２６】
続いて、採用する統計的推定法について説明する。データ分布の統計モデルに含まれるパラメータθの確率分布におけるハイパーパラメータαや、状態量ｚの確率分布に内在するハイパーパラメータｗ（各不良要因ｘ１…の確率分布パラメータや応答曲面モデルに内在するパラメータなど）を考慮する場合には、ハイパーパラメータα、ｗ の事前分布を決定する（ステップ５５、５６）。例えばデータ分布のハイパーパラメータの事前分布の決定の場合には、観測データ集合が従う確率分布へ事前に当てはめた結果から得られた分布パラメータの期待値と信頼区間をもとに事前分布を決定することが考えられる。
【００２７】
ハイパーパラメータの事前分布の考慮により、経験に基づく知識などの主観的情報を反映することが可能となる。ここで、考慮するのはα、ｗのどちらかのみ、あるいはどちらも考慮しない場合もある。ハイパーパラメータを考慮する場合は、データ分布の統計モデルにおけるパラメータθの推定は確率分布として得られ、考慮しない場合は点推定となる。ハイパーパラメータの推定には、ＡＢＩＣやＥＩＣなどの情報量基準に基づくＭＡＰ推定などの統計的推定手法を用いることができる。
【００２８】
次に、各不良モードの各仮説について、状態量ｚ１…の確率分布をそれぞれ求める（ステップ５７、ステップ３１〜３９、材料／形状／負荷データベース４１、ＣＡＥ解析モデルベース、応答曲面モデルベース４３、観測データベース４４）。状態量の確率変数集合は、対象とする数値実験点について、ＣＡＥによる現象解析の結果得られた状態量ｚ１…のデータ集合をもとに、不良要因ｘ１…と状態量ｚ１…の因果相関関係を応答曲面モデルにより表現した上で、応答曲面モデルを活用したモンテカルロシミュレーションにより算出する。ここで、連続値の状態量については、分割数と各分割範囲を指定した上で離散化を行い、状態量の確率分布の離散化を行った上で、確率密度が０でないデータ集合を抽出する。これにより、階層型統計モデル推定時の計算数を大幅に合理化できる
【００２９】
次に、対象とする観測データ集合と、上記で設定したデータ分布、パラメータθの確率分布、状態量ｚの確率分布をもとに、データ分布の統計モデルに含まれるパラメータθを推定する（ステップ５８）。パラメータθの推定には、最尤法やＥＭ（Expectation Maximization）法を活用する。
【００３０】
得られたパラメータθの推定結果をもとに、データ分布と観測データの適合度評価（誤差評価）を行う（ステップ５９）。適合度評価スペックを満たすケースがただ一つの仮説のみの場合は仮説が同定できたことになる。適合度評価スペックを満たすケースが存在しない場合には、対象とする観測データの決定やデータ分布の決定やハイパーパラメータの事前分布の決定や状態量の確率分布の算出を修正して再度繰り返す。これらの繰返しでも適合する仮説が無い場合には、仮説の想定を修正して再度実施する。
【００３１】
次に、得られたデータ分布と状態量ｚ１…の確率分布から、ベイズ法により事後分布を算出する（ステップ６０）。事後分布の評価（ステップ６１）により、推定スペックを満たさない場合には、対象とする観測データの決定やデータ分布の決定やハイパーパラメータの事前分布の決定や状態量の確率分布の算出を修正して再度繰り返す。これらの繰返しでも適合する仮説が無い場合には、仮説の想定を修正して再度実施する。
【００３２】
次に、推定結果の安定性を評価する（ステップ６２）。上記の推定結果は、変数の離散化条件やサンプリングの選び方などの推定条件に影響を大きく受ける可能性があり、その推定条件の変動に対する安定性を評価しておくことは重要である。推定条件を変更した上で、再度推定を行うことを繰返し、推定条件の変動に対する安定性（ロバスト性）を評価することができる。安定性が確認できれば、各仮説に対する状態量データを決定する（ステップ６３）。
【００３３】
適合する仮説が唯一得られた場合には、適合仮説について得られた状態量データ集合を満たす不良要因ｘの組合せを、応答曲面モデルを用いて探索することにより、不良要因データの組み合わせを同定する。適合する仮説が複数得られた場合は、以下の不良要因データの同定を行う。以下では状態推定法と多重仮説検定を用いた場合を示す。
【００３４】
まず、不良要因データの推定のための応答曲面モデルを状態空間表現する（ステップ８１）。この状態空間表現をもとに、状態方程式と観測方程式を表現する（ステップ８２）。上記工程にて得られた対象とする状態量データ集合を正規変数化し、各仮説に対して状態量の観測データとする（ステップ８３）。
【００３５】
次に、対象とする不良要因ｘ１…が従う確率分布から各仮説に対する不良要因のサンプリングデータを生成する（ステップ８４）。各仮説に対応する応答曲面モデルから対象とする状態量データのサンプリングデータを不良要因のサンプリングデータについてモンテカルロ法により生成する（ステップ８５）。
【００３６】
カルマンフィルタなどの状態推定法により、各仮説に対応する観測データをもとに、不良要因データの推定と逐次更新および残差分布の推定を行う（ステップ８６）。
【００３７】
各仮説に対応する応答曲面モデルの多重仮説検定により、最も適合するモデルを抽出し、不良発生メカニズム（適合仮説）を同定する（ステップ８７）。得られた応答曲面モデルから、不良要因データの組合せを同定する（ステップ８８）。ここで、同定した仮説に対して、不良発生確率が小さくなるような不良要因の組合せを探索してもよい（ステップ８９）。
【００３８】
以上説明のように、材料／製造／試験／負荷／品質モニタリングによるデータや現象解析結果が蓄積されるに従い、不良現象に関するモデル情報や不良メカニズム情報の未知の部分に対して、積極的に新たな因果相関構造を有する統計モデルで表現することは、有用と考えられる。このようなモデルにより、能動的に情報抽出や不良メカニズムの早期検出・同定を行うための方法論を構築することができる。
【００３９】
図４は、以上説明した半導体装置設計支援方法を実施するための装置の概略構成を示す図である。図４において、図３中に示したものと対応するものには同一符号を付してある。この装置１００は、ＣＡＥによる前処理部１０１、統計モデルの推定部１０２、不良発生仮説の同定部１０３を内部に有し、さらに、入力装置１１０、出力装置１２０、ＲＯＭ１３０、ＲＡＭ１４０、材料／形状／負荷データベース４１、ＣＡＥ解析モデルベース４２、応答曲面モデルベース４３、観測データベース４４、データ分布・故障物理モデルベース７１を周辺に有する。
【００４０】
ＣＡＥによる前処理部１０１は図３のステップ３１〜３９を実行し、統計モデルの推定部１０２は図３のステップ５１〜６３を実行し、不良発生仮説の同定部１０３は図３のステップ８１〜８９を実行する。入力装置１１０は、装置１００を操作するためのユーザ入力を与えるマンマシンインターフェースである。出力装置１２０は、ユーザにこの装置１００に関する各種の情報を表示、出力するためのマンマシンインターフェースである。ＲＯＭ１３０およびＲＡＭ１４０は、この装置１００を制御するためのソフトウエアの格納や、処理で必要な記憶領域の提供などを行うものである。
【００４１】
以下では、より具体的な例として、半導体パッケージの熱変形メカニズム（不良メカニズム）の同定問題を対象として、半導体装置設計支援方法の具体的な流れを説明する。
【００４２】
素子の薄型化および材料と形状に内在する相対的なばらつきの増大に起因して、半導体パッケージの熱変形はＭ字、Ｗ字、鞍型、凸型など、材料特性の組み合わせのわずかな差により複雑な変形挙動を示す。材料特性や寸法のばらつきだけではなく、測定データにも大きなばらつきが含まれており、熱変形メカニズムを明確化するには、不良現象に関する解析・評価技術の改良だけではなく、ＣＡＥ結果と実測データを活用して、統計的に同定するための新たな手法が必要である。
【００４３】
試験・検査などの実験データやシミュレーションによる数値実験データなどの現実の観測モニタリングデータに対して情報処理を行う場合、実験計画やノイズ除去を十分行った上で観測を行った場合を除いては、はじめから少数のパラメータで統計モデルを構成できるケースは少ない。不確定要素をすべて除去できるように試験・観測法や評価技術を改良することが重要であるが、開発期間や費用対効果の面や技術的な問題でパラメータの少数化が困難な場合も多い。本例では、複数の不良現象に関する観測データを統計的手法により情報処理する際に、潜在変数を含む多変量データを、複合的な現象を説明するのに必要な多くのパラメータをすべて含んだ形で推定することに主眼を置く。この種類の問題を取り扱う方法としてベイズ法がある。まず、ベイズ法による情報処理の構成要素に即して全体の枠組みを説明する。
【００４４】
観測データとデータ全体の確率分布から尤度関数を構成できる場合、統計モデルのパラメータを推定するのに最尤法がある。しかしモデルのパラメータの数が多く、観測データ数が不十分な場合にはうまく推定できない場合もある。このような場合、従来は、多項式回帰モデルなどを活用して“Reparametrization”を実施しパラメータを少数化してきたが、得られた多項式のパラメータに実用的な意味を見出せないことが多い。
【００４５】
ベイズ法のポイントは、２つの確率分布の同時分布が表現できれば、その一方の確率変数が観測された時点で、もう一方の確率分布が推定できるというものである。同時分布を構成するにあたって、データ分布と事前分布を適切に想定することが重要となる。本例では、従来、合理的な想定が困難な場合が多かった事前分布を、ＣＡＥと応答曲面法およびモンテカルロ法を用いたシーケンスにより各不良要因の不確定性に対する複数の潜在変数に関する確率分布を算出することにより設定した。また、データ分布は、不良データの発生機構に関する知識を含んだ上で想定する必要がある。本例では、不良現象の故障物理モデルを考慮した信頼度基底関数を含んだ形で、各潜在変数の影響度の把握が容易な比例ハザードモデルをそれぞれの不良モードに対して想定した。
【００４６】
［対象］
対象とする半導体パッケージのＣＡＥ解析モデルを図５に示す。この解析モデルにより熱変形解析を行う。対象とする不良モードは、１）半導体パッケージのコプラナリティ（反り量）と、２）はんだ接合部の熱疲労寿命、である。また、対象とする観測データは、１）半導体パッケージの反り分布データ（ｘ方向とｙ方向の反り分布を解析データとする）、２）温度サイクル試験の熱疲労寿命データ、である（図９、図１０、図１１を参照）。
【００４７】
状態量としては以下のように、上記の不良現象において重要な物理量である１２変数を対象とした（図８を参照）。反り分布測定データのみでは、データ量不足とノイズに起因して、半導体パッケージの熱変形の曲率や水平方向の変位といった状態量の値を抽出することはできないため、これらの変数は、直接には観測できない潜在変数（状態量）と言える。
【００４８】
図６は、半導体パッケージの熱変形に寄与すると考えられる要因を示す図である。熱変形要因は複雑に絡み合っており、同じ変形モードを表現するための仮説が複数考えられる。反りメカニズムを明確化するには、シミュレーションと実測モニタリング（およびフィルタリング）を組み合わせた統計的な同定方法が不可欠である。図７は、半導体パッケージの１／４の部分における熱変形モードの例を示す図である。剥離の存在位置や部材の材料特性の組合せなどで反りモードは大きく変化する。
【００４９】
図１２は、想定した仮説（熱変形メカニズム）を示す図である。また、図１３は、本例で対象とする統計モデルを示す図である。対象とする不良データｙとして、１）半導体パッケージの反り分布、２）はんだ接合部の熱疲労寿命分布を設定し、さらに、不良要因ｕを含む形で統計モデル（階層ベイズモデル、応答曲面モデル）を構成した。状態量としては、すでに述べたように熱変形の挙動において重要な物理量である１２変数を対象とした（図８、図１５を参照）。
【００５０】
本例では、観測データに関するデータ分布として、熱疲労寿命データに対しては、ワイブル分布の信頼度基底関数をＲ0とした比例ハザードモデルを設定した。反り分布（ｘ方向、ｙ方向)は、今回の対象の場合、状態量が決まると反り分布を一意に規定できるとしているため、データ分布としては、反り量の測定点に関するデータ分布としてデルタ関数を設定した。このデータ分布の統計モデル（確率モデル）には、観測変数ｙ（今回の場合熱疲労寿命データ）、潜在変数（状態量）ｚおよびモデルパラメータθ（ワイブル分布パラメータρ、βおよび比例ハザードモデルパラメータのＢ）が含まれる。
【００５１】
データ分布は以下のように表される（ただし状態量の数が１２のためｋ＝１２）
【数８】
2007108843.tif 000010
ここで、Ｂは統計モデルのパラメータであり、観測データをもとに推定される。また、Ｒ0とＲ（ｙ｜ｚ）は以下の式から算出できる。
【数９】
2007108843.tif 000011
【数１０】
2007108843.tif 000012
【数１１】
2007108843.tif 000013
この比例ハザードモデルには、Ｒ0に含まれる未知パラメータと、Ｂiの未知パラメータが内在する。
【００５２】
ワイブル分布の信頼度関数は以下のとおり表される。
【数１２】
2007108843.tif 000014
ここで、データ分布の統計モデルに含まれるパラメータθの確率分布におけるハイパーパラメータαの事前分布は、観測データ集合が従う確率分布へ事前に当てはめた結果から得られた分布パラメータの期待値と信頼区間をもとに事前分布を決定し、情報量基準に基づくＭＡＰ推定によりパラメータ値を決定した。
【００５３】
状態量の確率分布：
次に、各不良モードの各仮説について、状態量の確率分布をそれぞれ求める。状態量の確率変数集合は、対象とする数値実験点について、ＣＡＥによる現象解析の結果得られた状態量のデータ集合をもとに、不良要因と状態量の因果相関関係を応答曲面モデルにより表現した上で、応答曲面モデルを活用したモンテカルロシミュレーションにより算出する。
【００５４】
ここで、連続値の状態量については、分割数と各分割範囲を指定した上で離散化を行い、状態量の確率分布の離散化を行った上で、確率密度が０でないデータ集合を抽出した。今回の場合、不良要因ｕについて正規分布を設定し、平均値±２０％の範囲でモンテカルロシミュレーションを行い、各状態量について最大値と最小値の間を１０分割することにより離散化した上で、確率密度が０でない状態量のデータ集合を算出した。大規模なモンテカルロシミュレーションとなるため、より細分割化して高精度化するには、分散並列化などによる計算処理の高速化の必要性があると考えられる。
【００５５】
応答曲面モデルのパラメータ値：
構造解析ソフトＡＢＡＱＵＳＶVer6.4を用いた熱変形シミュレーション結果を利用して、実験計画法Ｌ１８に基づく応答曲面法により、１２個の状態量についてそれぞれ応答曲面モデルを作成した。１２個の応答曲面モデルに含まれるパラメータ値を仮説２の場合について図１４に示す。ここでｕは等価な不良要因を示している。状態量ｚは、図１５に示す通りである。また、参考のためＬ１８の実験計画法に基づく分散分析による各要因の影響度解析結果を図１６、図１７に示す（仮説２の場合）。
【００５６】
統計モデルパラメータの推定：
対象とする観測データ集合と、上記で設定したデータ分布、パラメータθの確率分布、状態量ｚの確率変数データ集合をもとに、データ分布の統計モデルに含まれるパラメータθ（今回の対象の場合、ワイブル分布パラメータρ、βおよび比例ハザードモデルパラメータのＢ）を各仮説について推定する。今回のパラメータθの推定にはＥＭ法を活用した。
【００５７】
統計的同定による仮説検定：
本例では、カルマンフィルタと多重仮説検定を用いた仮説検証を試みた（図１８、図１９、図２０、図２１、図２２を参照）。まず、各仮説の応答曲面モデルを状態空間表現した（図１９を参照）。
【００５８】
同定した仮説に対応する統計モデルの推定にて得られたパラメータθをもとに、データ分布と観測データの適合度評価（誤差評価）を行った結果、最大誤差は±１０％以下であった。今回対象とした半導体パッケージの場合、熱変形測定結果と信頼性試験結果を再現するには、インターポーザ表層の材料特性の異方性と、パッケージ樹脂の応力緩和による見かけの弾性率の低下および等価線膨張率の増加を考慮する必要があることが明らかになった。
【００５９】
以下では、本実施形態で用いるパラメータ推定法などの統計的手法について補足説明する。
【００６０】
［ＭＡＰ推定］
ｙというデータに基づいて推定した将来観測値zの分布をＰ（ｚ｜ｙ）で表す。Ｐ（ｚ｜ｙ）の良さをその平均対数尤度∫Ｑ（ｚ）ｌｏｇＰ（ｚ｜ｙ）ｄｚの大きさで評価することができる。Ｑ（ｚ）は確率変数Ｚの一般に未知の真の分布である。平均対数尤度が大きいモデルほど真の分布Ｑ（ｚ）に対する近似が良いと考えられる。この基準を情報量基準という。情報量基準の観点から事後分布で平均化した分布の最適性を比較できる。データｙとｚの発生機構Ｑ（ｙ）がＰ（ｙ｜θ＊）という関数で書ける場合を考える。以下ではθの分布π（θ＊）がわかっている場合について考える。
【００６１】
情報量基準ＡＢＩＣ：
赤池らは以下の情報量基準ＡＢＩＣを定義しハイパーパラメータθをＡＢＩＣの最小化によって決定する方法を提案した。第２項目はθの次元による補正項であり、データ分布や事前分布のモデル比較・選択を可能にするためのものである。
すなわち、情報量基準ＡＢＩＣ（α，ｗ）は、
【数１３】
2007108843.tif 000015
で表される。
【００６２】
Ｐ（ｚ｜ｙ）の平均対数尤度は
【数１４】
2007108843.tif 000016
と表現できる。
特定のｚの場合のみに良い推定値となる推定法では困るため、上記の平均対数尤度をデータｚの分布で平均化すると以下で表される。
【数１５】
2007108843.tif 000017
θそのものは未知のため、θの分布で平均化すると以下の式になる。
【数１６】
2007108843.tif 000018
【００６３】
最後の式が最大値をとるのは以下のときである。
【数１７】
2007108843.tif 000019
対数関数の形から任意のｓ＞０に対してｌｏｇ（ｓ）≦ ｓ−１
が成立、これから
【数１８】
2007108843.tif 000020
【数１９】
2007108843.tif 000021
【数２０】
2007108843.tif 000022
【００６４】
パラメータの真の事前分布がわかっている場合、ベイズ法で求められる事後分布のもとでの混合分布：
【数２１】
2007108843.tif 000023
がデータｙが観測された場合の最適な予測分布となる。π（θ｜ｙ）単峰でとがっていればそれを最大とするθの値：
【数２２】
2007108843.tif 000024
をパラメータの値として持つ以下の
【数２３】
2007108843.tif 000025
が近似的な予測分布となり、パラメータの推値としては
【数２４】
2007108843.tif 000026
がＭＡＰ推定値となる。
【００６５】
パラメータの真の事前分布が既知であれば、情報量基準の観点からＭＡＰ推定の利用が正当化される。パラメータの真の事前分布が既知というのは厳しい条件であるが、仮定したパラメトリック事前分布の族｛π（θ｜α，ｗ）｝が真の事前分布を含んでいる場合にはＡＢＩＣ最小化法で選んだ事前分布モデルに基づくＭＡＰ推定が正当化されると考えてよい。
【００６６】
情報量基準ＥＩＣ：
｛π（θ｜α，ｗ）｝が真の事前分布を含んでいる場合にＭＡＰ推定が情報量基準の意味でよい推定を与える。また、｛π（θ｜α，ｗ）｝が真の事前分布を含んでいない場合でも、ＭＡＰ推定モデルの平均対数尤度：
【数２５】
2007108843.tif 000027
が他の選択肢、例えば非ベイズ的なモデルの対数尤度：
【数２６】
2007108843.tif 000028
より良かった場合には、その推定値を使うのが良いと考えられる。
【００６７】
情報量基準ＡＩＣは以下の式から得られる。すなわち、ＡＩＣは
【数２７】
2007108843.tif 000029
である。第２項のＳとＴはｙとは独立であるが、ｙと同一の分布に従う確率変数である。この部分が
【数２８】
2007108843.tif 000030
がデータｙにもとづくθの最尤推定値であろうとなかろうと第２項はｙと同じ分布に従う互いに独立な乱数の集合｛Ｓｉ｝と｛Ｔｉ｝を用いて
【数２９】
2007108843.tif 000031
と近似計算できる。
【００６８】
｛Ｓｉ｝と｛Ｔｉ｝としてブートストラップサンプルｙ＊を利用するのがＥＩＣである。
ＥＩＣは以下の式で定義される情報量基準である。すなわち、ＥＩＣは
【数３０】
2007108843.tif 000032
である。この場合、
【数３１】
2007108843.tif 000033
が最尤推定量でなければならないということはなく、ＭＡＰ推定でよい。ＡＩＣの比較でモデル選択ができるのと同じように、ＥＩＣの比較によって同じモデルパラメータの選定法が選択できることになる。
【００６９】
ＥＩＣの実際の計算は次に式による。
【数３２】
2007108843.tif 000034
ただし、ｙＹｉ＝｛ｙＹ１，ｙＹ２，…，ｙＹＫ｝がブートストラップサンプルである。ＥＩＣの計算にあたっては、生データのリサンプリングではなく、残差のリサンプリングのほうが有用であることもある。
【００７０】
［ハイパーパラメータの選択と識別問題］
ハイパーパラメータの値が異なれば、当然のことながら同じ観測データから異なるモデルが推定される。事前情報をできる限り活用するために事前分布を導入するが、事前分布がハイパーパラメータを含む場合にはその選定法を与える必要がある。この問題を系統的に解決するために赤池らは最尤法に基づく情報量基準を活用する方法を提案している。
【００７１】
潜在変数の事前分布とデータ分布から、データと潜在変数の同時分布が決まる。データ分布や事前分布がハイパーパラメータ（本例の場合はデータ分布に含まれるパラメータαや不良要因の想定確率分布および応答曲面モデルのパラメータｗ）に依存する場合には、同時分布がハイパーパラメータに依存することになる。形式的にデータと潜在変数の組をデータの分布のモデルと考えれば、データを観測し潜在変数を欠損値とした上で、最尤法でハイパーパラメータを推定することができる。
【００７２】
潜在変数（欠損値）に対しては、潜在変数の周辺を積分し周辺分布として潜在変数ｚをintegrate-outするか、もしくは仮の値で代用する。ハイパーパラメータの推定においては、ＡＢＩＣやＥＣＩなどの情報量基準に基づくＭＡＰ推定が利用できる。データ分布のパラメータθの推定には、計算効率に優れたＥＭ推定を用いることが可能である。本例では、ｗは確定値として扱っており、αのみをハイパーパラメータとしてＭＡＰ推定法により推定した上で、ＥＭ法によりパラメータθを推定している。ｗもハイパーパラメータとして考慮する場合には、各ｗの更新値ごとに事前分布を算出し直した上で情報量基準を求める必要がある。また、ブートストラップサンプリングを行う場合は、観測データｙのみ対象とする場合と、観測データｙと潜在変数ｚの両方についてサンプリングを繰り返す場合の両方が可能である。ブートストラップ法に基づく情報量基準ＥＩＣは、事前分布モデルが真の分布を含まない場合のＭＡＰ推定の利用を可能にする。
【００７３】
［最尤推定］
統計モデルは対象とする現象を表現する数理モデルの一つであり、現象について仮定した統計モデルの未知パラメータを推定し、データの平滑化、同定、予測を行うことが、統計的推定の目的である。データを規定する変数に不規則性があり、その不確定性をある確率分布に従う確率変数として表現する点が統計モデルの特徴である。
【００７４】
統計モデルの変数は、直接データとして観測される観測変数と本来直接観測不可能な潜在変数とに大別される。潜在変数は、状態量、画素、欠測値、隠れ変数とも呼ばれることがある。潜在変数を導入することによって、より柔軟なモデルを構成することができる。潜在変数を含む統計モデルは数多く存在する。例えば、混合分布推定問題の場合、観測データはどの要素分布から生成されたかは本来観測できない。他にも音声認識などで用いられる隠れマルコフモデル、時系列解析で用いられる状態空間モデルなども潜在変数をもつモデルである。ＥＭ法は、これら潜在変数をもつ統計モデルのモデルパラメータの最尤推定値を求めるための一般的数値解法である。Dempter、Laird、Rubinによって１９７７年に考案されたものであり、グラフィカルモデリングや共分散構造分析などのパラメータ推定などにも活用されている。
【００７５】
不完全データの最尤推定：
観測データをＤ、観測データのそれぞれに対応する潜在変数からなる集合をＺとする。このとき、Ｄを不完全データ、（Ｄ，Ｚ）を完全データとよぶ。潜在変数が離散値をとるものとし、完全データの確率分布がｐ（Ｄ，Ｚ；θ）と与えられている場合を考える。ただしθは未知パラメータを表している。観測データ（不完全データ）に対する対数尤度関数は以下の式で表現できる（潜在変数が連続値の場合は、和を積分に置き換えればよい）。
【数３３】
2007108843.tif 000035
【００７６】
Ｚに関して周辺化されているので、潜在変数を含む統計モデルにおける上記尤度関数は周辺尤度関数ともよばれる。観測データＤが与えられたもとでの未知パラメータθの最尤推定値は対数尤度関数Ｌを最大化するパラメータ値：
【数３４】
2007108843.tif 000036
として定義される。ここでargmax(*)は*を最大化する変数を表す。最尤推定値とは直観的には、観測データＤを最も支持するパラメータ値である。ただし、最尤推定法ではその支持度を対数尤度関数で評価する。上記最大化は非線形最適化問題に属し解析的に解くことができないため、ニュートン法などの非線形最適化手法を利用する必要がある。ＥＭ法はより効率的に推定するための解析法である。
【００７７】
［ＥＭ法］
ＥＭ法では上記の対数尤度関数を直接最大化する代わりに、より簡単な式になる場合の多い完全データの対数尤度関数の条件付き期待値の最大化を行う。ＥＭ法は２種類のステップの逐次反復により最尤推定値を求める。第ｔ反復でのパラメータの推定値をθ（ｔ）とすると、まず期待値計算ステップ（Ｅステップ）で完全データ対数尤度関数：
【数３５】
2007108843.tif 000037
の条件付き期待値：
【数３６】
2007108843.tif 000038
を計算する。ここで記号Ｅ｛ｆ（ｘ｜ｙ）｝はｙが与えられたときのｆ（ｘ）の条件付き期待値を表す。
【数３７】
2007108843.tif 000039
はＤを観測したもとでのＺの事後分布で、ベイズの定理より次式で計算される。
【数３８】
2007108843.tif 000040
【００７８】
ＥテップでのＱ関数の計算は、直感的には、現在のパラメータ推定値に基づいて算出される潜在変数の事後分布：
【数３９】
2007108843.tif 000041
をその信頼度として不完全情報を暫定的に補っていると解釈できる。次に、最大化ステップ（Ｍステップ）でＱ関数を最大化するθを第ｔ＋１反復でのパラメータの推定値とする。これらＥＭステップを収束条件が満たされるまで繰り返す。
【００７９】
ＥＭ法ではＱ関数を計算する際、すべての可能な潜在変数の値に関する和を計算する。このとき、潜在変数の取り得る値の数が指数オーダになると、その和を現実的な時間で計算することが困難となる。そこで本例では、潜在変数ｚを不良要因に関する応答曲面モデルで表現し、応答曲面モデルに基づくモンテカルロシミュレーションにより潜在変数に関する確率分布において、確率密度が０とならない潜在変数の組合せに関するデータセットを事前に抽出することにより、計算規模の縮小化を行う。
【００８０】
観測変数ｙ＝｛ｙ１，ｙ２，…，ｙＭ｝、潜在変数Ｚ＝｛ｚ１，ｚ２，…，ｚＭ｝とするとき、観測データに対する尤度関数は以下のように書ける。ここで、θ＝（θy，θz）はモデルパラメータを表している。
【数４０】
2007108843.tif 000042
ＥステップにおけるＱ関数は、
【数４１】
2007108843.tif 000043
と書ける。
（ｙ，Ｚ）の同時分布は次式のようになる。
【数４２】
2007108843.tif 000044
【００８１】
また、ベイズの定理から次式が計算できる。
【数４３】
2007108843.tif 000045
このようなＥステップの計算に対して、Ｅステップをサンプリングにより近似するモンテカルロ法が適用できる。
【００８２】
ベイズ推定は最尤推定と異なり、予測値ではなく予測分布を求める。また、ベイズ推定は、観測データが少ない場合、最尤推定に対し汎化能力の高い推定が可能となる利点をもつ。このように、ベイズ推定では、データを観測したもとでの事後分布が中心的な役割を果たす。ベイズ推定ではすべての未知量が確率変数として取り扱われる。未知パラメータθ、モデル構造Ｍ、 潜在変数集合Ｚをもつ統計モデル（確率モデル）に対しては、通常、以下に示すグラフィカルモデルが仮定される。この場合、観測データＤはＺが与えられたもとでは、モデル構造Ｍに独立（ｐ（Ｄ｜Ｚ，Ｍ）＝ｐ（Ｄ｜Ｚ））であることを示す。
【数４４】
2007108843.tif 000046
【００８３】
上記より、同時分布は次式にように分解できる。
【数４５】
2007108843.tif 000047
ここで、ｐ（Ｄ，Ｚ｜θ，Ｍ）はモデル構造が既知のもとでの完全データ（Ｄ，Ｚ）の尤度関数を表し、ｐ（θ｜Ｍ）、ｐ（Ｍ）はパラメータ、モデル構造の事前分布をそれぞれ表す。上式より観測データ集合Ｄを得たもとでのすべての未知量の同時事後分布はベイズの定理より以下の式で与えられる。
【数４６】
2007108843.tif 000048
【００８４】
上式の分母はｐ（Ｄ）と等価である。またパラメータθの事後分布はＺ、Ｍに関する周辺化より
【数４７】
2007108843.tif 000049
として得られる。事後分布が得られれば、新たなサンプルｄに対する事後予測分布は次式により得ることができる。
【数４８】
2007108843.tif 000050
上式は変数を離散化することによりモンテカルロシミュレーションを活用して近似的に求めることができる。
【図面の簡単な説明】
【００８５】
【図１】複数の設計変数を入力変数とし複数の不良モードにかかわる測定可能量を出力変数とし、かつ、それらの間に複数の潜在的な状態変数が介在するとした階層構造の統計モデルを説明するための概観図。
【図２】一般的な故障物理モデルを例示する説明図。
【図３】本発明の一実施形態に係る半導体装置設計支援方法を示すフロー図。
【図４】図３に示した半導体装置設計支援方法を実施するための装置の概略構成を示す図。
【図５】半導体パッケージの熱変形メカニズム（不良メカニズム）の同定問題における対象とする半導体パッケージのＣＡＥ解析モデルを示す図。
【図６】半導体パッケージの熱変形メカニズム（不良メカニズム）の同定問題における不良要因を示す図。
【図７】半導体パッケージの熱変形メカニズム（不良メカニズム）の同定問題における半導体パッケージの１／４の部分における熱変形モードの例を示す図。
【図８】半導体パッケージの熱変形メカニズム（不良メカニズム）の同定問題における状態量を示す図。
【図９】半導体パッケージの熱変形メカニズム（不良メカニズム）の同定問題における観測データのひとつを示す図。
【図１０】半導体パッケージの熱変形メカニズム（不良メカニズム）の同定問題における観測データの別のひとつを示す図。
【図１１】半導体パッケージの熱変形メカニズム（不良メカニズム）の同定問題における観測データのさらに別のひとつを示す図。
【図１２】半導体パッケージの熱変形メカニズム（不良メカニズム）の同定問題における想定した仮説（熱変形メカニズム）を示す図。
【図１３】半導体パッケージの熱変形メカニズム（不良メカニズム）の同定問題における対象とする統計モデルを示す図。
【図１４】半導体パッケージの熱変形メカニズム（不良メカニズム）の同定問題において求められた応答曲面モデルの例を示す図。
【図１５】図１３に示した統計モデルにおける状態量を示す図。
【図１６】半導体パッケージの熱変形メカニズム（不良メカニズム）の同定問題において求められた各要因の影響度解析結果を示す図。
【図１７】図１６の続図であって、半導体パッケージの熱変形メカニズム（不良メカニズム）の同定問題において求められた各要因の影響度解析結果を示す図。
【図１８】半導体パッケージの熱変形メカニズム（不良メカニズム）の同定問題における、カルマンフィルタと多重仮説検定を用いた仮説検証の流れを示す図。
【図１９】各仮説の応答曲面モデルの状態空間表現を示す図。
【図２０】等価な不良要因の探索範囲を示す図。
【図２１】半導体パッケージの熱変形メカニズム（不良メカニズム）の同定問題における、カルマンフィルタと多重仮説検定による仮説検証を示す図。
【図２２】半導体パッケージの熱変形メカニズム（不良メカニズム）の同定問題における同定結果を示す図。
【符号の説明】
【００８６】
１１…入力変数、１２…状態変数、１３…出力変数、１４…ＣＡＥデータ、１５…観測データ、１６…故障物理モデル、４１…材料／形状／負荷データベース、４２…ＣＡＥ解析モデルベース、４３…応答曲面モデルベース、４４…観測データベース、７１…データ分布・故障物理モデルベース、１００…半導体装置設計支援装置、１０１…ＣＡＥによる前処理部、１０２…統計モデルの推定部、１０３…不良発生仮説の同定部、１１０…入力装置、１２０…出力装置、１３０…ＲＯＭ、１４０…ＲＡＭ。
【図１】
2007108843.tif 000051
【図２】
2007108843.tif 000052
【図３】
2007108843.tif 000053
【図４】
2007108843.tif 000054
【図６】
2007108843.tif 000055
【図７】
2007108843.tif 000056
【図８】
2007108843.tif 000057
【図９】
2007108843.tif 000058
【図１０】
2007108843.tif 000059
【図１１】
2007108843.tif 000060
【図１２】
2007108843.tif 000061
【図１３】
2007108843.tif 000062
【図１４】
2007108843.tif 000063
【図１５】
2007108843.tif 000064
【図１６】
2007108843.tif 000065
【図１７】
2007108843.tif 000066
【図１８】
2007108843.tif 000067
【図１９】
2007108843.tif 000068
【図２０】
2007108843.tif 000069
【図２１】
2007108843.tif 000070
【図２２】
2007108843.tif 000071
【図５】
2007108843.tif 000072