声門逆フィルタ法

声門逆フィルタ法^[1]（せいもんぎゃくフィルタほう、英: glottal inverse filtering）は、有声音の音圧信号または口腔気流信号から声道の影響を推定して取り除き、声門を通過する体積流量の波形、またはその時間微分に対応する波形を推定する音声分析法である^[2]^:623–625^[3]^:53。音圧信号を入力とする場合には、声道に加えて口唇放射の影響を補償する^[2]^{:623, 625–626}^[3]^:53–54。日本語文献では「声門逆フィルタリング」という表現も用いられる^[4]^:642。推定対象は文献で声門流、声門励振または声門音源と呼ばれるが、推定される量が体積流量に対応する波形か、その時間微分に対応する波形かは、入力信号と補償処理によって異なる^[2]^:624–625^[3]^{:53, 57}。

声門逆フィルタ法は単一のアルゴリズムではなく、声道特性を分析者が調整する方法、声門の閉鎖期を利用する方法、線形予測法を反復的または重み付きで用いる方法、二次計画法を用いる方法、音源と声道を同時に推定する方法、信号の位相特性を利用する方法、深層学習を用いる方法などを含む方法群である^[5]^{(著者最終稿PDFの15–21頁)}^[3]^:57。推定波形とそこから算出される特徴量は、音声生成と声質の研究、音声合成、話者認識、感情認識、音声に基づく疾患分類の研究などに用いられている^[2]^:637–640^[5]^{(著者最終稿PDFの21–27頁)}^[3]^:61–67。ただし、得られる波形は声門流の直接測定値ではなく、音声生成モデル、録音系、分析条件およびアルゴリズムの仮定に依存する推定値である^[2]^:640–641。

原理

有声音の生成は、第一近似として、声帯振動によって生じる声門音源が声道で共鳴し、口唇から放射される過程として表される^[2]^:624–626^[3]^:53–54。線形のソース・フィルタモデルでは、観測音声と声門音源のZ変換をそれぞれ $S(z)$ 、 $G(z)$ 、声道と口唇放射の伝達関数をそれぞれ $V(z)$ 、 $L(z)$ とすると、 $S(z)=G(z)V(z)L(z)$ と表す^[2]^:625–626^[3]^:53–54。声門逆フィルタ法では、声道と口唇放射の伝達特性を表すモデルの逆特性を観測信号へ適用し、 $G(z)$ の推定値を求める^[2]^:625–626。

マイクロホンで口唇前方の音圧を記録する場合、口唇放射は口唇における体積流量を放射音圧へ変換する微分器に近い特性を持つ^[2]^:624–626^[3]^:53–54。このため、声道の逆フィルタリングに加えて口唇放射を補償し、実装によっては逆フィルタ後の信号を積分して声門流に対応する波形を得る^[2]^{:625–626, 640}。マイクロホンや録音装置の低周波域における振幅特性または位相特性が不十分であると、積分後の波形形状に誤差が生じ得る^[2]^:640。

別の入力として、口唇周囲に装着したマスクによって口腔気流を記録する方法がある^[2]^:627。Martin Rothenbergが1973年に報告した方式では、周囲に通気部を設けた流量測定用マスクを使用する^[6]^:1632。この方式では、逆フィルタ後の信号を絶対流量に校正でき、直流成分を含む声門流を扱える一方、装置の装着を要し、測定可能な帯域にも制約がある^[2]^:627。電気声門図は声帯接触に関連する信号を与え、声門閉鎖時点の推定を補助できるが、声門流そのものを測定するものではない^[2]^:628^[5]^{(著者最終稿PDFの5, 15–16頁)}。

推定には、音源と声道に関する仮定が必要となる^[2]^{:625, 628–630}^[3]^:53–57。多くの方法は、声道を全極型フィルタで近似すること、声門閉鎖期では音源の寄与が小さいこと、声門音源と声道で位相特性が異なること、または声門流を所定のモデルで表せることなどを利用する^[2]^:628–630^[5]^{(著者最終稿PDFの15–18頁)}^[3]^:55–57。実際の発声では声門音源と声道の相互作用が生じるため、音源とフィルタを独立した線形系として分けるモデルは近似である^[3]^:55–56^[7]^:98。

歴史

Alkuによる2011年の総説は、録音音声から声道の影響を打ち消して声門音源を推定した初期の研究として、Robert L. Millerが1959年に発表した研究を挙げている^[8]^:667^[2]^:626。Millerの装置はアナログ回路で第1フォルマントの影響を補償し、声門閉鎖後の推定流量がゼロとなるように逆フィルタを調整するものであった^[2]^:626。1960年代にはスウェーデン王立工科大学を中心に、複数のフォルマントを扱うアナログ逆フィルタが研究され、閉鎖期の波形に残るフォルマント振動を最小化することが調整基準として用いられた^[2]^:627。

1968年、Alan V. OppenheimとRonald W. Schaferは、畳み込みで結合した音源と声道の成分をケプストラム領域で加算成分へ変換する準同型解析を音声へ適用した^[9]^:221–222^[2]^:627。1970年、M. NakatsuiとJ. Suzukiは、時間領域のディジタル逆フィルタによって声門音源波形を観測する方法を報告した^[10]^:664–665^[2]^:627。1970年代には、口腔気流を入力とするRothenbergの方法と、閉鎖期の標本から線形予測によって声道を推定する閉鎖期共分散法が報告された^[6]^:1632^[11]^:350^[2]^:627–628。

1980年代には、声門音源の入力モデルと声道の自己回帰モデルを同時に推定する方法が研究された^[2]^:629。Paul Milenkovicは1986年、線形の音源入力モデルと自己回帰系を同時に推定する方法を報告した^[12]^:28^[2]^:629。1992年にはPaavo Alkuが、低次の線形予測による声門成分の推定と高次の線形予測による声道推定を段階的に行う反復適応逆フィルタリングを発表した^[13]^:109^[2]^:630。

2000年代以降には、逆フィルタリングと声門流モデルへの適合を同時に行う方法、Z変換の零点配置または複素ケプストラムを用いて因果成分と反因果成分を分ける方法、重み付き線形予測によって声道を推定する方法などが提案された^[14]^:479^[15]^:344–345^[16]^:855–857^[17]^:596。2010年代後半以降には、既存の逆フィルタ法または音声合成器から得た教師信号を用いて、音声から声門流を推定する深層学習型の方法も報告された^[18]^:95–96^[19]^:1327–1328。

推定法

手動調整

手動調整では、分析者が音声スペクトルを参照して各フォルマントに対応する逆フィルタの中心周波数と帯域幅を設定し、推定波形とそのスペクトルを確認しながら調整する^[2]^:626–627。調整基準の一つは、声門が閉じていると考えられる区間に声道共鳴による振動が残らず、波形が平坦または滑らかになることである^[2]^:627–628。個々の信号に合わせて設定できる一方、中心周波数と帯域幅の試行的な調整には時間を要し、設定は分析者の判断に依存する^[2]^:626–628^[3]^:56。

閉鎖期を用いる方法

閉鎖期を用いる方法は、声門が閉じていると推定される区間では声門音源の寄与が小さく、観測信号が主として声道の自由応答を表すという仮定を用いる^[2]^:628。閉鎖期共分散法では、選択した区間へ共分散基準による線形予測を適用し、全極型の声道モデルを求める^[11]^:350^[2]^:628。閉鎖期の位置は音声信号から推定するほか、電気声門図から得た声門閉鎖時点を補助的に利用できる^[2]^:628。

この方法は、閉鎖期が明確で十分な長さを持つ音声では声道推定に利用できるが、閉鎖期の開始・終了位置の誤差に敏感である^[2]^:628。基本周波数が高い音声または息もれ声では利用できる閉鎖期が短くなり、線形予測に必要な標本数を確保しにくい^[2]^{:628, 641}。

二次計画法を用いる方法

Airaksinenらは、閉鎖期分析の原理を利用しながら、声道と口唇放射を一つのフィルタで表し、その係数を二次計画法によって求めるQPR法を提案した^[20]^:929。QPR法の最適化では、推定波形の閉鎖期を平坦にする条件を課すことができ、原著では合成された持続母音と実音声を用いて評価された^[20]^:929。

反復適応逆フィルタリング

反復適応逆フィルタリング（iterative adaptive inverse filtering、IAIF）は、低次の線形予測で声門音源の大まかなスペクトル傾斜を推定し、その影響を軽減した信号から高次の線形予測で声道を推定する処理を段階的に行う^[13]^:109^[2]^:630。個々のフォルマントを手動調整せずに処理できる一方、声道モデルの線形予測次数、声門流モデルの線形予測次数および口唇放射係数の設定が推定結果へ影響する^[21]^:24–27。

IAIFには、一階線形予測で一つの実極をモデル化し、その影響を逆フィルタリングで除去する処理を、自己回帰係数の絶対値が所定の閾値を下回るまで反復するIOP-IAIFと、入力信号を積分した後に一階線形予測と逆フィルタリングを3回行い、声門流モデルの線形予測次数を3次に固定するGFM-IAIFなどの派生形がある^[21]^:26^[22]^:4–5。Mokhtariらの比較では、標準IAIF、IOP-IAIF、GFM-IAIFの誤差は、母音、基本周波数、発声条件およびパラメータ設定によって異なった^[21]^{:24, 28–36}。

音源と声道の同時推定

音源と声道を同時に推定する方法では、声門音源に所定の入力モデルを置き、声道モデルとともにパラメータを求める^[2]^:629–630^[3]^:57。Milenkovicの方法は、自己回帰系として表した声道と線形入力モデルを同時に推定する^[12]^:28。SIMは、逆フィルタリングと声門流モデルへの適合を同時に行い、声道フィルタと声門流モデルのパラメータを求める^[14]^:479。

SahooとRoutrayは2016年、連結管で表した声道とLiljencrants–Fantモデルで表した声門流の時間微分を状態空間モデルへ組み込み、拡張カルマンフィルタによって連続音声から状態とモデルパラメータを推定する方法を報告した^[23]^:1230。AlzamendiとSchlotthauerは2017年、状態空間モデルのパラメータ推定に期待値最大化法を用い、声門音源と声道フィルタを同時に推定する方法を報告した^[24]^:5^[3]^:57。

位相特性を用いる方法

位相特性を用いる方法は、有限長の音声区間を、主として最小位相となる成分と、反因果または最大位相の成分へ分ける混合位相モデルを利用する^[15]^:344–345^[16]^:855–857。Z変換の零点を用いる方法では、零点の単位円に対する位置から因果成分と反因果成分を分離する^[15]^:344–346。複素ケプストラムを用いる方法では、複素ケプストラムの正および負のケフレンシーを利用して分解する^[16]^:855–859。

これらの混合位相分解法は、声門開放時点および声門閉鎖時点の同定と、分析窓の長さおよび窓関数の選択に影響される^[19]^:1327。Drugmanらの比較研究は、閉鎖期法、IAIFおよび混合位相分解を合成音声と連続発話で比較し、清浄条件では閉鎖期法または混合位相法が良好である一方、雑音条件ではIAIFの相対的な頑健性が高かったと報告した^[25]^:20。

重み付き線形予測

QCP法で用いられる重み付き線形予測では、線形予測誤差に標本ごとの重みを与え、声門音源の主励振が声道モデルの推定へ及ぼす影響を抑える^[17]^:596。quasi closed phase（QCP）法は、複数の基本周期にわたる音声標本を対象とし、声門閉鎖時点に基づくattenuated main excitation（AME）重み関数によって閉鎖期の標本を相対的に重視して声道モデルを推定する^[17]^:596^[22]^{:2, 5}。QCP法は信号全体の標本を使用して閉鎖期を重視するため、単一の短い閉鎖区間だけを抽出する方式ではない^[22]^:2。AME重み関数の形状は、声門閉鎖時点とposition quotient（PQ）、duration quotient（DQ）、ramp quotient（RQ）によって定められる^[22]^:5。重み付き線形予測には、QCPのほかにも重み関数または誤差基準が異なる方法が提案されている^[3]^:57。

深層学習を用いた推定

深層学習を用いた推定では、音声信号から声門流または声門流の時間微分を推定するようにニューラルネットワークを学習する^[18]^:95^[19]^:1327–1328。自然音声について基準となる声門流を大規模に得ることが難しいため、既存の逆フィルタ法で得た推定値、または音声合成器で生成した音声と既知の声門流の組を教師信号として用いる研究がある^[18]^:95^[19]^:1327–1328。

Narendraらは2019年、符号化電話音声から抽出した特徴量を、対応する非符号化音声にQCP法を適用して得た声門流へ写像する深層ニューラルネットワークを報告した^[18]^:95。Langheinrichらは2022年、調音合成器VocalTractLabによって生成した連続音声と声門流の組を用いて双方向長短期記憶ネットワークを学習し、音声から声門流の時間微分を推定する方法を報告した^[19]^:1327–1329。同研究は、自然音声に対応する基準波形を必要な規模で得られないことを、合成データを用いる理由としている^[19]^:1327–1328。

推定波形のパラメータ化

声門逆フィルタ法で得た波形は、そのまま波形形状を観察するほか、一周期の時間構造、振幅およびスペクトル形状を特徴量で表すためにパラメータ化される^[2]^:631–635^[5]^{(著者最終稿PDFの18–21頁)}。時間構造を表す特徴には、基本周期に占める声門開放期間の割合である声門開放率（open quotient、OQ）、声門流が増加する区間と減少する区間の長さの比であるspeed quotient（SQ）、基本周期に占める減少区間の割合であるclosing quotient（ClQ）などがある^[2]^:631–632。これらは開放時点、閉鎖時点および最大流量時点の定義に依存し、推定波形に声道共鳴の残留振動または雑音があると算出が難しくなる^[2]^:631–633。

振幅に基づく特徴には、声門流の交流振幅、最小流量、声門流の時間微分の負の最大振幅（maximum flow declination rate、MFDR）などがある^[2]^:632–634。絶対流量と直流成分を評価するには、流量マスクなどによる校正済みの測定が必要となる^[2]^{:627, 640}。normalized amplitude quotient（NAQ）は、声門流の交流振幅を、基本周期と声門流の時間微分の負の最大振幅の積で正規化した量である^[26]^:701–703^[2]^:632。

周波数領域の特徴には、第1調波と第2調波の振幅差（H1–H2）、スペクトル傾斜、基本波の振幅に対する上位調波の振幅和の比であるharmonic richness factor（HRF）などがある^[2]^:634–635^[5]^{(著者最終稿PDFの21頁)}。parabolic spectral parameter（PSP）は、推定声門流のスペクトル減衰を放物線近似によって定量化する指標として提案された^[27]^:67^[2]^:635。文献上は、推定声門流のスペクトルから求めるH1–H2と、放射音声に声道補正を施して求める同種の指標があり、算出対象となる信号と補正方法を区別する必要がある^[5]^{(著者最終稿PDFの21頁)}。

推定波形をLiljencrants–Fantモデルなどのパラメトリックな声門音源モデルへ適合させ、モデルパラメータで波形を表す方法もある^[5]^{(著者最終稿PDFの5–6, 20頁)}。声門音源モデルはモデルごとにパラメータ数と表現する波形特性が異なるため、得られるパラメータの意味は使用したモデルに依存する^[5]^{(著者最終稿PDFの5–6, 20頁)}。

評価方法

自然発話では、分析対象の音声に対応する声門体積流量の基準波形を非侵襲的に得ることが難しい^[2]^:640–641^[28]^:38–41。このため、評価には既知の声門音源と声道フィルタから生成した合成音声、声帯振動と声道音響を計算する物理モデル、生理学的または調音的な音声合成器、人工音源で励振する3Dプリント声道模型、マイクロホン信号・電気声門図・声帯高速度映像を同時記録した自然音声などが用いられる^[2]^:640–641^[29]^:1718^[21]^:24–28^[28]^:38–43。

既知の音源波形と声道フィルタから音声を生成する評価では、生成に用いた音源波形を基準として使用できる^[2]^:640–641^[28]^:38–41。一方、音声生成モデルと評価対象の逆フィルタ法が同じ線形ソース・フィルタモデルなどの仮定を共有する場合には、評価対象に有利な条件となり得る^[2]^:640–641^[21]^:25^[28]^:39–41。物理モデルや生理学的な音声合成器は、同じ単純な線形ソース・フィルタ仮定だけに依存した評価を避けるために用いられる^[2]^:641^[21]^:25^[28]^:39–41。一方、これらの信号も声帯振動、声道形状または音響伝搬のモデルから生成されたものであり、自然発話に対応する直接測定値ではない^[29]^:1718^[21]^:25–28^[28]^:38–43。OPENGLOTは、複数の生成方式による評価用データを共通環境で提供するものとして2019年に発表された^[28]^:38–43。同環境の自然音声データには基準となる声門流は含まれず、音圧信号、電気声門図および声帯高速度映像が収録されている^[28]^{:38, 42}。

評価指標には、基準波形と推定波形の波形誤差のほか、OQ、closing quotient（ClQ）、NAQ、quasi-open quotient（QOQ）、H1–H2、HRF、スペクトル傾斜などの特徴量誤差が用いられる^[28]^:41–42^[22]^{:3, 7–10}。一つの比較でも複数の波形指標と特徴量指標を併用することがある^[29]^:1718^[22]^:7–10。雑音を加えた条件、基本周波数、母音、発声様式または分析パラメータを変えた条件での頑健性も評価される^[25]^:20^[21]^{:24, 28–36}^[22]^{:1, 7, 9–13}。

Drugmanらは2012年、閉鎖期法、IAIFおよび混合位相分解を合成音声と連続発話で比較し、清浄条件と雑音条件で方式の相対的な結果が変化したと報告した^[25]^:20。Chienらは2017年、生理学的な調音音声合成器による持続母音と連続発話を用いて複数方式を比較し、声質、基本周波数、声門下圧、母音および発話形式に応じた誤差を評価した^[29]^:1718。Mokhtariらは2018年、計算物理モデルを用いて三つのIAIF系方式を比較した^[21]^:24。Freixesらは2023年、OPENGLOTで生成した男性音声と女性音声を模した合成母音を用いてIAIF系とQCP系の方式を比較した^[22]^:1。

応用

発声と声質の研究

声門逆フィルタ法は、声門流波形と音声の音響特性を対応付け、発声様式、声質、感情、声の強さ、話者差および歌声などを調べるために用いられている^[2]^:637–640^[5]^{(著者最終稿PDFの21–24頁)}。推定波形のOQ、NAQ、スペクトル傾斜などを用いて、息もれ声、通常発声、圧迫発声などの発声様式に伴う声門音源の違いが研究されている^[2]^:637–639^[5]^{(著者最終稿PDFの21–24頁)}。これらの特徴量の算出値は、逆フィルタ法、録音条件、基本周波数および音圧などの影響を受ける^[2]^{:631–635, 637–641}^[5]^{(著者最終稿PDFの18–24頁)}。Kreimanらは、声道補正後のH1*–H2*、OQおよび声門面積波形の非対称性の関係が話者間で大きく異なったと報告している^[30]^:2625。

流量マスクを用いる研究では、声門流の直流成分、交流振幅、最小流量および声門流の時間微分の負の最大振幅などを絶対値として扱えるため、声の強さなどに伴う気流変化の研究に利用されてきた^[2]^{:627, 632–634, 637–640}。マイクロホンで記録した音声から得た推定値は、流量校正を伴わない限り絶対流量または直流流量としては扱えず、主として波形形状または正規化特徴量の分析に用いられる^[2]^:640。

音声技術

音声技術では、推定した声門音源と声道フィルタを分けてモデル化し、音声合成、声質変換および音声符号化へ用いる研究がある^[5]^{(著者最終稿PDFの24–27頁)}。Raitioらは2011年、声門逆フィルタ法で分離した音源と声道の特徴をそれぞれ統計的にモデル化する隠れマルコフモデル（HMM）音声合成方式を報告した^[31]^:153。推定した声門特徴を話者認識、感情認識および発話様式の分析へ加える研究も報告されている^[2]^:637–640^[5]^{(著者最終稿PDFの21–27頁)}。

疾患分類の研究

音声から疾患を自動分類する研究では、声門逆フィルタ法から得た特徴量を統計解析または自動分類の入力として用いる例がある^[2]^:639–640^[3]^:61–67。2026年の総説章は、声門特徴のみを用いる研究、メル周波数ケプストラム係数など他の音響特徴と組み合わせる研究、推定した声門流を深層学習分類器へ直接入力する研究を扱っている^[3]^:61–67。

同章は、深層学習による分類には多くの学習データが必要である一方、疾患を対象とする音声データセットは小規模であることが多いと述べている^[3]^:67–68。また、音声に基づく健康状態のバイオマーカーは臨床診断に代わるものではないと明記している^[3]^:59–60。

限界と誤差要因

自然発話では基準となる声門流を非侵襲的に得ることが難しいため、個々の推定波形を既知の基準波形と直接比較することは困難である^[2]^:640–641^[28]^:38–41。音源と声道を独立した線形系として扱う逆フィルタ法では、両者の相互作用に由来する非線形性はモデルに含まれない^[3]^:55–56。PalaparthiとTitzeの計算モデルでは、音源–声道相互作用の増加に伴って平均二乗誤差は低下した^[7]^:98。一方、基本周波数が第1フォルマントへ近い条件では、声道共鳴周波数の推定が難しくなり、予測誤差が増加した^[7]^:98。

基本周波数が高い音声では調波間隔が広くなり、声道のスペクトル包絡を推定する情報が減少する^[2]^:641^[3]^:55–56。閉鎖期法では一周期と閉鎖期が短くなり、線形予測に利用できる標本数も減る^[2]^{:628, 641}。息もれ声または不完全閉鎖を伴う発声では、明瞭で十分な長さの閉鎖期を前提とする方法の適用が難しくなる^[2]^:628^[29]^:1718。

鼻音および鼻音化母音では、多くの方式が用いる全極型声道モデルだけでは声道の周波数特性を十分に表せず、推定精度が低下し得る^[2]^:641。連続発話では、声道形状、基本周波数および発声様式が時間的に変化し、持続母音より分析条件が複雑になる^[2]^:641–642^[29]^:1718。

背景雑音は推定結果に影響し、方式によって清浄条件と雑音条件で相対的な結果が変わる^[25]^:20。マイクロホンで記録した音声を積分して声門流に対応する波形を求める場合には、録音系の低周波域における振幅特性と位相特性も波形形状に影響し得る^[2]^:640。手動方式では分析者によるフォルマント周波数と帯域幅の調整が必要となる^[2]^:626–628。自動方式でも、IAIFでは声道・声門流モデルの線形予測次数と口唇放射係数、QCPでは声門閉鎖時点とPQ・DQ・RQ、混合位相分解では分析窓の長さと窓関数など、方式ごとの設定または推定が必要となる^[21]^:24–28^[22]^:5^[19]^:1327。比較研究では、同じ方式でも設定値によって誤差が変化することが報告されている^[21]^:28–36^[22]^{:4–5, 11–13}。