numpy.linalg.linalg.LinAlgError: Singular matrix のエラー

#問題点
Pythonを用いて、こちらのサイトを参考にしながら類似楽曲検索システムを制作しています。
制作していると、numpy.linalg.linalg.LinAlgError: Singular matrixというエラーが発生してしまいました。
どうすれば正常に動作するでしょうか？

#エラーメッセージ

bash
1python mir.py ../testdata/sig/Recollections.sig ../testdata/sig/ ../testdata/html/
2Traceback (most recent call last):
3  File "mir.py", line 175, in <module>
4    emd = calcEMD(targetSigPath, sigPath)
5  File "mir.py", line 68, in calcEMD
6    dist[i * numFeatures + j] = symKLDiv(mu1, S1, mu2, S2)
7  File "mir.py", line 42, in symKLDiv
8    return 0.5 * (KLDiv(mu1, S1, mu2, S2) + KLDiv(mu2, S2, mu1, S1))
9  File "mir.py", line 30, in KLDiv
10    invS2 = np.linalg.inv(S2)
11  File "/home/kei/.local/lib/python2.7/site-packages/numpy/linalg/linalg.py", line 528, in inv
12    ainv = _umath_linalg.inv(a, signature=signature, extobj=extobj)
13  File "/home/kei/.local/lib/python2.7/site-packages/numpy/linalg/linalg.py", line 89, in _raise_linalgerror_singular
14    raise LinAlgError("Singular matrix")
15numpy.linalg.linalg.LinAlgError: Singular matrix
16

#コードの内容

python
1#coding:utf-8
2import os
3import sys
4import numpy as np
5import numpy.linalg
6import rpy2.robjects as robjects
7from collections import defaultdict
8
9# mir.py
10# usage: python mir.py [sig file] [sig dir] [html file]
11# sig file  : クエリ楽曲のシグネチャファイル
12# sig dir   : 検索対象のシグネチャファイルのディレクトリ
13# html file : 検索結果を出力するHTMLファイル
14
15# 引数で指定したシグネチャファイルに近い
16# 上位N件の楽曲を出力する
17
18
19import numpy as np
20import numpy.linalg
21
22def KLDiv(mu1, S1, mu2, S2):
23    """正規分布間のカルバック・ライブラー情報量"""
24    # 逆行列を計算
25    try:
26        invS1 = np.linalg.inv(S1)
27    except numpy.linalg.linalg.LinAlgError:
28        raise;
29    try:
30        invS2 = np.linalg.inv(S2)
31    except numpy.linalg.linalg.LinAlgError:
32        raise;
33
34    # KL Divergenceを計算
35    t1 = np.sum(np.diag(np.dot(invS2, S1)))
36    t2 = (mu2 - mu1).transpose()
37    t3 = mu2 - mu1
38    return t1 + np.dot(np.dot(t2, invS2), t3)
39
40def symKLDiv(mu1, S1, mu2, S2):
41    """対称性のあるカルバック・ライブラー情報量"""
42    return 0.5 * (KLDiv(mu1, S1, mu2, S2) + KLDiv(mu2, S2, mu1, S1))
43    
44import rpy2.robjects as robjects
45
46# Rで輸送問題を解くライブラリ
47# Rのデフォルトパッケージではないのでインストールが必要
48# Rでinstall.packages("lpSolve")
49robjects.r['library']('lpSolve')
50transport = robjects.r['lp.transport']
51
52def calcEMD(sigFile1, sigFile2):
53    # シグネチャをロード
54    sig1 = loadSignature(sigFile1)
55    sig2 = loadSignature(sigFile2)
56
57    # 距離行列を計算
58    numFeatures = sig1.shape[0]                 # クラスタの数
59    dist = np.zeros(numFeatures * numFeatures)  # 距離行列（フラット形式）
60
61    for i in range(numFeatures):
62        mu1 = sig1[i, 1:21].reshape(20, 1)   # 縦ベクトル
63        S1 = sig1[i, 21:421].reshape(20, 20)
64        for j in range(numFeatures):
65            mu2 = sig2[j, 1:21].reshape(20, 1)
66            S2 = sig2[j, 21:421].reshape(20, 20)
67            # 特徴量iと特徴量j間のKLダイバージェンスを計算
68            dist[i * numFeatures + j] = symKLDiv(mu1, S1, mu2, S2)
69
70    # シグネチャの重み（0列目）を取得
71    w1 = sig1[:,0]
72    w2 = sig2[:,0]
73
74    # 重みと距離行列からEMDを計算
75    # transport()の引数を用意
76    costs = robjects.r['matrix'](robjects.FloatVector(dist),
77                                 nrow=len(w1), ncol=len(w2),
78                                 byrow=True)
79    row_signs = ["<"] * len(w1)
80    row_rhs = robjects.FloatVector(w1)
81    col_signs = [">"] * len(w2)
82    col_rhs = robjects.FloatVector(w2)
83    
84    t = transport(costs, "min", row_signs, row_rhs, col_signs, col_rhs)
85    flow = t.rx2('solution')
86    
87    dist = dist.reshape(len(w1), len(w2))
88    flow = np.array(flow)
89    work = np.sum(flow * dist)
90    emd = work / np.sum(flow)
91    return emd
92
93
94
95def loadSignature(sigFile):
96    """シグネチャファイルをロード"""
97    mat = []
98    fp = open(sigFile, "r")
99    for line in fp:
100        line = line.rstrip()
101        mat.append([float(x) for x in line.split()])
102    fp.close()
103    return np.array(mat)
104
105def getArtist(mp3Path):
106    """MP3ファイルからアーティストを取得"""
107    import eyeD3
108    try:
109        tag = eyeD3.Tag()
110        tag.link(mp3Path)
111        artist = tag.getArtist()
112    except:
113        artist = "None"
114    # 空白のとき
115    if artist == "": artist = "None"
116    return artist
117
118def makeHTML(ranking, htmlFile, N=10):
119    """ランキングをHTML形式で出力"""
120    import codecs
121    fout = codecs.open(htmlFile, "w", "utf-8")
122
123    # HTMLヘッダを出力
124    fout.write('<!DOCTYPE html>\n')
125    fout.write('<html lang="ja">\n')
126    fout.write('<head><meta charset="UTF-8" /><title>%s</title></head>\n' % htmlFile)
127    fout.write('<body>\n')
128    fout.write('<table border="1">\n')
129    fout.write(u'<thead><tr><th>ランク</th><th>EMD</th><th>タイトル</th>')
130    fout.write(u'<th>アーティスト</th><th>音声</th></tr></thead>\n')
131    fout.write(u'<tbody>\n')
132
133    # ランキングを出力
134    rank = 1
135    for sigFile, emd in sorted(ranking.items(), key=lambda x:x[1], reverse=False)[:N]:
136        prefix = sigFile.replace(".sig", "")
137
138        # rawをwavに変換（HTMLプレーヤー用）
139        rawPath = os.path.join("raw", prefix + ".raw")
140        wavPath = os.path.join("wav", prefix + ".wav")
141        if not os.path.exists("wav"): os.mkdir("wav")
142        os.system('sox -r 16000 -e signed-integer -b 16 "%s" "%s"' % (rawPath, wavPath))
143
144        # アーティスト名を取得
145        mp3Path = os.path.join("mp3", prefix + ".mp3")
146        artist = getArtist(mp3Path)
147
148        # HTML出力
149        # HTML5のオーディオプレーヤーを埋め込む
150        audio = '<audio src="%s" controls>' % wavPath
151        fout.write("<tr><td>%d</td><td>%.2f</td><td>%s</td><td>%s</td><td>%s</td></tr>\n"
152                   % (rank, emd, prefix, artist, audio))
153        rank += 1
154
155    fout.write("</tbody>\n");
156    fout.write("</table>\n")
157    fout.write("</body>\n")
158    fout.write("</html>\n")
159    fout.close()
160
161if __name__ == "__main__":
162    if len(sys.argv) != 4:
163        print "python mir.py [sig file] [sig dir] [html file]"
164        sys.exit()
165
166    targetSigPath = sys.argv[1]
167    sigDir = sys.argv[2]
168    htmlFile = sys.argv[3]
169
170    ranking = defaultdict(float)
171    
172    # 全楽曲との間で距離を求める
173    for sigFile in os.listdir(sigDir):
174        sigPath = os.path.join(sigDir, sigFile)
175        emd = calcEMD(targetSigPath, sigPath)
176        if emd < 0: continue
177        ranking[sigFile] = emd
178
179    # ランキングをEMDの降順にソートして出力
180    N = 10
181    rank = 1
182    for sigFile, emd in sorted(ranking.items(), key=lambda x:x[1], reverse=False)[:N]:
183        print "%d\t%.2f\t%s" % (rank, emd, sigFile)
184        rank += 1
185
186    # EMDの昇順に上位10件をHTMLにして出力
187    makeHTML(ranking, htmlFile, N)
188