[PYTHON] Verschieben von Objekttrennungen mit Robust PCA
Einführung
Im vorherigen Artikel (https://qiita.com/matsxxx/items/652e58f77558faecfd23) habe ich ein Beispiel für das Trennen von sich bewegenden Objekten in einem Video mit Dynamic Mode Decomposition (DMD) gezeigt. In diesem Artikel werden wir versuchen, die Objekttrennung mit der Robust Principal Component Analysis (Robust PCA) des PCP-Algorithmus (Principal Component Pursuit) zu verschieben [1].
Umgebung
- Windows 10 home edition
- Python (3.7.6)
- Numpy(1.18.5)
- Opencv(4.1.1)
Robust PCA via PCP
Robustes PCA ist eine Methode zur Zerlegung in eine Struktur mit niedrigem Rang und eine Struktur mit geringer Dichte. Im Fall eines sich bewegenden Körpers hat der Hintergrund eine niedrige Rangstruktur und der sich bewegende Körper eine spärliche Struktur. Dieses Mal werde ich Robust PCA ausprobieren, einen PCP-Algorithmus. Zerlegen Sie die ursprüngliche Datenmatrix $ M $ in eine niedrigrangige Matrix $ L $ und eine spärliche Matrix $ S $.
M = L + S
PCP zerlegt sich in $ L $ und $ S $, indem der folgende erweiterte Lagrange $ l $ minimiert wird.
l(L, S, Y) = \| L \|_* + \lambda\|S\|_1 +<Y, M - L - S> + \frac{\mu}{2}\|M-L-s\|^2_F
$ L $ ist eine niedrigrangige Matrix, $ S $ ist eine dünn besetzte Matrix, $ Y $ ist eine Lagrange-Multiplikatormatrix und $ \ lambda $ und $ \ mu $ sind PCP-Parameter. $ \ <, * > $ Ist das innere Produkt. Die Dimensionen jeder Variablen sind $ M, L, S, Y \ in \ mathbb {R} ^ {n_1 \ times n_2}
\begin{align}
&||A||_* = \mathrm{tr}(\sqrt{A^\mathrm{T}A})= \sum^{\mathrm{min}(n_1,n_2)}_{i=1} \sigma_i\\
&||A||_1 = \underset{1\leq j \leq n_2}{\mathrm{max}}\sum^{n_1}_{i=1}|a_{ij}|\\
&||A||_F = \sqrt{\sum^{n_1}_{i=1}\sum^{n_2}_{j=1}|a_{ij}|^2} = \sqrt{\mathrm{tr}(A^\mathrm{T}A)} = \sqrt{\sum^{\mathrm{min}(n_1,n_2)}_{i=1}\sigma^2_i}
\end{align}
Es wird sein. Der Index $ \ mathrm {T} $ repräsentiert eine transponierte Matrix. Die Formel, die die Norm beschreibt, lautet [Wikipedia-Beschreibung](https://ja.wikipedia.org/wiki/%E8%A1%8C%E5%88%97%E3%83%8E%E3%83%AB% Ich bezog mich auf E3% 83% A0). Die Sparse-Matrix $ S $ wird mit dem Schrumpfungsoperator aktualisiert, und die niedrigrangige Matrix $ L $ wird mit dem Singularwert-Schwellenwertoperator aktualisiert. Der Schrumpfungsoperator $ S_ \ tau $ des Parameters $ \ tau $ ist
\begin{align}
&S_\tau(x) = \mathrm{sgn}(x) \mathrm{max}(|x|-\tau, 0)
\end{align}
Es wird ausgedrückt als. Der Singularwert-Schwellenwertoperator $ D_ \ tau $ verwendet den Schrumpfungsoperator $ S_ \ tau $ und die Singularwertzerlegung.
D_{\tau}(X) = US_\tau(\Sigma)V^\mathrm{T}, \
\mathrm{where} \ X = U\Sigma V^\mathrm{T}
Es wird sein. $ U, \ Sigma und V $ sind der linke Singularvektor, der Singularwert und der rechte Singularvektor von $ X $. Verwenden Sie den Schrumpfungsoperator $ S_ \ tau $ und den Singularwert-Schwellenwertoperator $ D_ \ tau $, um die niedrigrangige Matrix $ L $ und die dünn besetzte Matrix $ S $ zu aktualisieren. Die Aktualisierungsformel wird als $ X_k $ vor der Aktualisierung der Variablen $ X $ und $ X_ {k + 1} $ nach der Aktualisierung geschrieben, und der PCP-Parameter $ \ tau, \ mu $ wird eingegeben und geschrieben.
\begin{align}
&L_{k+1} = D_{1/\mu}(M - S_k - \mu^{-1}Y_k) \\
&S_{k+1} = S_{\lambda/\mu}(M - L_{k+1} + \mu^{-1}Y_k)
\end{align}
Es wird sein. Das $ S $ -Symbol ist verwirrend, folgt jedoch der Notation in Referenz [1], mit der Ausnahme, dass die kontingente Matrix als transponierte Matrix umschrieben wird. Auch die Methode zum Einstellen der Parameter $ \ mu und \ lambda $ des Schrumpfungsoperators $ S_ \ tau $ und des Schwellenwertoperators $ D_ \ tau $ für Singularwerte ist in dem in arXiv platzierten PCP-Papier von 2009 falsch. , Seien Sie vorsichtig (2009 arXiv-Version ist $ D_ \ mu, S_ {\ lambda \ mu} $, aber es sollte $ D_ {1 / \ mu}, S_ {\ lambda / \ mu} sein Es sieht aus wie $). Die Lagrange-Multiplikatormatrix $ Y $ wird nach der erweiterten Lagrange-Methode berechnet.
Y_{k+1} = Y_k + \mu(M - L_{k+1} - S_{k+1})
Und aktualisieren. Die PCP-Parameter $ \ lambda $ und $ \ mu $ sind in Referenz [1] die Datenmatrix $ M \ in \ mathbb {R} ^ {n_1 \ times n_2} $.
\begin{align}
&\lambda = \frac{1}{\sqrt{\mathrm{max}(n_1,n_2)}} \\
&\mu = \frac{n_1 n_2}{4\|M\|_1}
\end{align}
Ist eingestellt. Das Folgende ist eine Zusammenfassung des PCP-Berechnungsverfahrens.
\begin{align}
&\mathrm{initialize:} S_0 = 0, Y_0 =0 \\
& \mathrm{while}\ \mathrm{not}\ \mathrm{converged}:\\
&\qquad L_{k+1} = D_{1/\mu}(M - S_k + \mu^{-1}Y_k) \\
&\qquad S_{k+1} = S_{\lambda/\mu}(M - L_{k+1} + \mu^{-1}Y_k) \\
&\qquad Y_{k+1} = Y_k + \mu(M - L_{k+1} - S_{k+1}) \\
&\mathrm{end} \ \mathrm{while} \\
&\mathrm{output:}L,S.
\end{align}
Die Konvergenz wird durch die folgende Frobenius-Norm bestimmt.
\|M -L - S\|_F \leq \delta\|M\|_F \quad \mathrm{with} \ \delta=10^{-7}
Code Robust PCA des oben genannten PCP-Algorithmus in Python.
import numpy as np
import numpy.linalg as LA
def rpca(M, max_iter=800,p_interval=50):
def shrinkage_operator(x, tau):
return np.sign(x) * np.maximum((np.abs(x) - tau), np.zeros_like(x))
def svd_thresholding_operator(X, tau):
U, S, Vh = LA.svd(X, full_matrices=False)
return U @ np.diag(shrinkage_operator(S, tau)) @ Vh
i = 0
S = np.zeros_like(M)
Y = np.zeros_like(M)
error = np.Inf
tol = 1e-7 * LA.norm(M, ord="fro")
mu = M.shape[0] * M.shape[1]/(4 * LA.norm(M, ord=1))
mu_inv = 1/mu
lam = 1/np.sqrt(np.max(M.shape))
while i < max_iter:
L = svd_thresholding_operator(M - S + mu_inv * Y, mu_inv)
S = shrinkage_operator(M - L + mu_inv * Y, lam * mu_inv)
Y = Y + mu * (M - L - S)
error = LA.norm(M - L - S, ord='fro')
if i % p_interval == 0:
print("step:{} error:{}".format(i, error))
if error <= tol:
print("converted! error:{}".format(error))
break
i+=1
else:
print("Not converged")
return L, S
Video
Vorheriger Artikel Verwenden Sie in ähnlicher Weise das Atrium des Datensatzes in hier. Ich werde. Ich benutze 120frame bis 170frame. Verwenden Sie Robust PCA, um die gehende Person vom Hintergrund zu trennen.
import cv2
def prepare_frames_to_use(video_path="./atrium_video.avi",#Bildpfad
start_frame=120,#Erster Frame
end_frame=170,#Letzter Frame
r=4,#Bildreduktionsparameter
):
#Bild wird geladen
cap = cv2.VideoCapture(video_path)
#Holen Sie sich Bildauflösung, Bildrate, Anzahl der Bilder
wid = cap.get(cv2.CAP_PROP_FRAME_WIDTH)#Seite
wid_r = int(wid/r)
hei = cap.get(cv2.CAP_PROP_FRAME_HEIGHT)#Vertikal
hei_r = int(hei/r)
fps = cap.get(cv2.CAP_PROP_FPS)#Bildrate
count = cap.get(cv2.CAP_PROP_FRAME_COUNT)#Anzahl der Frames
dt = 1/fps#Anzahl der Sekunden zwischen Frames
print(f"wid:{wid}", f" hei:{hei}", f" fps:{fps}", f" count:{count}", f"dt:{dt}")
#Extraktion des Zielrahmens
cat_frames = []
cap.set(cv2.CAP_PROP_POS_FRAMES,start_frame)
for i in range(end_frame - start_frame):
ret, img = cap.read()
if not ret:
print("no image")
break
buf = cv2.cvtColor(cv2.resize(img,(wid_r, hei_r)),
cv2.COLOR_BGR2GRAY).flatten()#
cat_frames.append(buf)
cap.release()
cat_frames = np.array(cat_frames).T
return cat_frames, wid_r, hei_r, fps
Moving Object Separation durch Robust PCA
Verwenden Sie Robust PCA, um Matrizen mit niedrigem Rang und geringer Dichte zu finden. Geben Sie einfach die Videodatenmatrix in die rpca-Funktion ein.
def detect_moving_object():
#Erfassung von zur Berechnung verwendeten Frames
frames, wid, hei, fps = prepare_frames_to_use()
#Führen Sie Robust pca aus
L, S = rpca(frames)
#Videoausgabe mit geringer Struktur
fourcc = cv2.VideoWriter_fourcc(*"mp4v")
writer = cv2.VideoWriter("sparse_video.mp4",
fourcc, fps, (int(wid), int(hei)))
for i in range(S.shape[1]):
s_buf = S[:,i]
#Helligkeitsanpassung
lum_min = np.abs(s_buf.min())
lum_max = np.abs(s_buf.max())
lum_w = lum_max + lum_min
s_buf = (s_buf + lum_min)/lum_w * 255
s_buf = cv2.cvtColor(s_buf.reshape(int(hei), -1).astype('uint8'),
cv2.COLOR_GRAY2BGR)
writer.write(s_buf)
writer.release()
Ergebnis
Wenn Sie sich das Video der spärlichen Prozession in der Abbildung rechts ansehen, können Sie sehen, dass die gehenden Menschen klar voneinander getrennt sind.
Verweise
- [1]: Candes, E. J., Li, X., Ma, Y., and Wright, J. 2011. Robust principal component analysis? J. ACM 58, 3, ` Article 11 (May 2011), 37 pages.
- [2]: Jodoin, J.-P., Bilodeau, G.-A., Saunier, N., Urban Tracker: Multiple Object Tracking in Urban Mixed Traffic, Accepted for IEEE Winter conference on Applications of Computer Vision (WACV14), Steamboat Springs, Colorado, USA, March 24-26, 2014.
Alle Quellcode
import numpy as np
import numpy.linalg as LA
def rpca(M, max_iter=800,p_interval=50):
def shrinkage_operator(x, tau):
return np.sign(x) * np.maximum((np.abs(x) - tau), np.zeros_like(x))
def svd_thresholding_operator(X, tau):
U, S, Vh = LA.svd(X, full_matrices=False)
return U @ np.diag(shrinkage_operator(S, tau)) @ Vh
i = 0
S = np.zeros_like(M)
Y = np.zeros_like(M)
error = np.Inf
tol = 1e-7 * LA.norm(M, ord="fro")
mu = M.shape[0] * M.shape[1]/(4 * LA.norm(M, ord=1))
mu_inv = 1/mu
lam = 1/np.sqrt(np.max(M.shape))
while i < max_iter:
L = svd_thresholding_operator(M - S + mu_inv * Y, mu_inv)
S = shrinkage_operator(M - L + mu_inv * Y, lam * mu_inv)
Y = Y + mu * (M - L - S)
error = LA.norm(M - L - S, ord='fro')
if i % p_interval == 0:
print("step:{} error:{}".format(i, error))
if error <= tol:
print("converted! error:{}".format(error))
break
i+=1
else:
print("Not converged")
return L, S
import cv2
def prepare_frames_to_use(video_path="./atrium_video.avi",#Bildpfad
start_frame=120,#Erster Frame
end_frame=170,#Letzter Frame
r=4,#Bildreduktionsparameter
):
#Bild wird geladen
cap = cv2.VideoCapture(video_path)
#Holen Sie sich Bildauflösung, Bildrate, Anzahl der Bilder
wid = cap.get(cv2.CAP_PROP_FRAME_WIDTH)#Seite
wid_r = int(wid/r)
hei = cap.get(cv2.CAP_PROP_FRAME_HEIGHT)#Vertikal
hei_r = int(hei/r)
fps = cap.get(cv2.CAP_PROP_FPS)#Bildrate
count = cap.get(cv2.CAP_PROP_FRAME_COUNT)#Anzahl der Frames
dt = 1/fps#Anzahl der Sekunden zwischen Frames
print(f"wid:{wid}", f" hei:{hei}", f" fps:{fps}", f" count:{count}", f"dt:{dt}")
#Extraktion des Zielrahmens
cat_frames = []
cap.set(cv2.CAP_PROP_POS_FRAMES,start_frame)
for i in range(end_frame - start_frame):
ret, img = cap.read()
if not ret:
print("no image")
break
buf = cv2.cvtColor(cv2.resize(img,(wid_r, hei_r)),
cv2.COLOR_BGR2GRAY).flatten()#
cat_frames.append(buf)
cap.release()
cat_frames = np.array(cat_frames).T
return cat_frames, wid_r, hei_r, fps
def detect_moving_object():
#Erfassung von zur Berechnung verwendeten Frames
frames, wid, hei, fps = prepare_frames_to_use()
#Führen Sie Robust pca aus
L, S = rpca(frames)
#Videoausgabe mit geringer Struktur
fourcc = cv2.VideoWriter_fourcc(*"mp4v")
writer = cv2.VideoWriter("sparse_video.mp4",
fourcc, fps, (int(wid), int(hei)))
for i in range(S.shape[1]):
s_buf = S[:,i]
#Helligkeitsanpassung
lum_min = np.abs(s_buf.min())
lum_max = np.abs(s_buf.max())
lum_w = lum_max + lum_min
s_buf = (s_buf + lum_min)/lum_w * 255
s_buf = cv2.cvtColor(s_buf.reshape(int(hei), -1).astype('uint8'),
cv2.COLOR_GRAY2BGR)
writer.write(s_buf)
writer.release()
if __name__ == "__main__":
detect_moving_object()