Alignment

Alignment

「對齊」。一堆數據,經過變換,盡量符合另一堆數據。已知甲堆數據、乙堆數據,求函數。

函數具有特殊限制,增進討論意義。以下討論相似變換。

Error

兩堆數據通常不一致。強硬地對齊,就會有「誤差」。

兩堆數據的誤差,有許多種衡量方式:

一、窮舉甲堆到乙堆的所有兩兩配對,累計距離。
二、窮舉甲堆每一筆數據,分別找到乙堆之中距離最近的那筆數據,累計距離。
三、窮舉乙堆每一筆數據,分別找到甲堆之中距離最近的那筆數據,累計距離。
四、二與三相加。
五、二與三聯集。
六、上述誤差,改為只累計前K短的距離。
七、上述誤差,只考慮互相鄰近的K筆數據。(部分銜接的效果)

隨隨便便就能發明一大堆誤差公式。重點在於計算時間長短、銜接密合程度。哪種誤差比較好?目前還沒有定論。

演算法(Principal Component Analysis)

兩堆數據各自計算PCA座標軸,然後對齊。

http://www.cs.tau.ac.il/~dcor/Graphics/cg-slides/svd_pca.pptx
http://www.cse.wustl.edu/~taoju/cse554/lectures/lect07_Alignment.pdf
http://perception.inrialpes.fr/~Horaud/Courses/pdf/Horaud_3DS_5.pdf

兩堆數據各自求得共變異矩陣,實施特徵分解。兩邊的特徵向量,依照特徵值大小一一對應。

最佳的位移量:數據中心的差距。最佳的旋轉量:特徵向量的夾角。最佳的縮放量:特徵值平方根的比例。

特徵向量不具方向性。指定每個特徵向量的方向,滿足右手座標系,才能求得旋轉矩陣,不含鏡射。

一種想法是窮舉各種可能的方向,從中找到旋轉角度最小者。

另一種想法是檢查determinant。det(E) = +1,是右手座標系。det(E) = -1,是左手座標系;只要再讓任意一個特徵向量顛倒方向,就是右手座標系。選擇最小的特徵值的特徵向量,讓誤差增加最少。

PCA座標軸容易因outlier而歪斜。是非常不精準的方法。

演算法(Procrustes Analysis)

當兩堆數據一樣多,並且預先知道數據對應關係,就有更精準的做法。

https://igl.ethz.ch/projects/ARAP/svd_rot.pdf
http://perception.inrialpes.fr/~Horaud/Courses/pdf/Horaud_3DS_6.pdf
    [ x1.x x2.x ... xN.x ]       [ y1.x y2.x ... yN.x ]
X = [ x1.y x2.y ... xN.y ]   Y = [ y1.y y2.y ... yN.y ]
    [ x1.z x2.z ... xN.z ]       [ y1.z y2.z ... yN.z ]
                                  _____________________ 
solve s R X + t = Y              | t: translate vector |
                            2    | s: scalar vector    |
minimize ‖ Y - (s R X + t) ‖     | R: rotation matrix  |
                                  ‾‾‾‾‾‾‾‾‾‾‾‾‾‾‾‾‾‾‾‾‾ 
then t = Ȳ - s R X̄
      2      T         T          2        2
     s = ∑ Ŷi Ŷi ÷ ∑ X̂i X̂i = ‖ Ŷ ‖F ÷ ‖ X̂ ‖F
            T                  T       T
     R = V U    with C₃ₓ₃ = X̂ Ŷ = U Σ V 

令函數是位移、旋轉、縮放三個步驟。方便起見,把位移改為最終步驟。數據套用函數之後,令誤差越小越好。誤差訂為對應數據的距離的平方的總和。

最佳的位移量:兩堆數據各自的平均值相減,相減之前先讓第一個平均值套用函數。最佳的縮放量:兩堆數據各自的變異數相除。最佳的旋轉量:兩堆數據各自減去平均值(中心位移至原點)、除以變異數(大小縮放為相等),求得維度的共變異矩陣,實施奇異值分解,然後VUᵀ矩陣相乘,即是旋轉矩陣。

R是正規正交矩陣,功效是旋轉暨鏡射。大家通常不想鏡射。若要避免鏡射,則必須滿足右手座標系。

det(R) = +1,是右手座標系。det(R) = -1,是左手座標系;只要再讓任意一個向量顛倒方向,就是右手座標系。選擇最小的特徵值所對應的向量,讓誤差增加最少。

演算法(Iterative Closest Point)

一、初始化:
 甲、實施PCA,求得對齊函數。
 乙、甲堆套用函數。讓甲堆靠近乙堆。
二、重複此步驟,直到最近點不變為止:
 甲、甲堆每一點各自找到在乙堆的最近點,建立對應關係。
   無視其餘點,實施PA,求得對齊函數。
 乙、甲堆套用函數。讓甲堆更加靠近乙堆。

尋找最近點,簡易的方式是窮舉法,進階的方式是乙堆套用分割區域的資料結構,諸如Uniform Grid、KD-Tree。請見本站文件「Region」。

演算法(RANSAC)

一、甲堆隨機抓幾點,乙堆隨機抓幾點,點數一樣多,推定它們依序一一對應。
  以這些點實施Procrustes Analysis,求得對齊函數。
二、計算「甲堆套用函數」與「乙堆」的誤差大小。
三、重複一二,取誤差最小者。

絕聖棄智,民利百倍。

Manifold Alignment

https://en.wikipedia.org/wiki/Manifold_alignment

Embedding(Under Construction!)

Embedding(Dimensionality Reduction)

「嵌入」。更換數據所在空間。例如從三維換成二維。例如從二維換成三維。例如從球面換成平面。例如從曲面換成環面。

嵌入方式有許多種。例如投影就是其中一種嵌入方式。例如訂立座標系統並且實施座標轉換,也是其中一種嵌入方式。

嵌入時,我們通常希望儘量保留數據的原始特性。特性可以是數據的聚散程度、數據之間的距離遠近。計算學家針對各種特性,設計各種演算法。

這裡先不談太廣,僅討論最簡單的情況:數據從N維空間換成M維空間。細分為三種情況:維度變小、維度不變、維度變大。

後兩者缺乏討論意義──數據不變、數據補零不就好了。因此我們僅討論維度變小。此時「嵌入」的意義就變成了:找到一個函數,讓一堆數據經過函數變換,減少數據維度。

演算法(Principal Component Analysis)

嵌入時,採用垂直投影。

實施Low-rank PCA,保留大的、捨棄小的特徵值暨特徵向量,只做投影、略過位移、略過縮放。

演算法(Metric Multidimensional Scaling)

嵌入時,盡量保留所有點對距離。

首先計算原數據X的兩兩距離,做為新數據Y的兩兩距離。而兩兩距離可以寫成「距離矩陣」:Mᵢⱼ = ‖pᵢ - pⱼ‖。

問題變成:已知距離矩陣M、求得數據Y。

Mᵢⱼ² = ‖qᵢ - qⱼ‖² = ‖qᵢ‖² + ‖qⱼ‖² - 2‖qᵢ ∙ qⱼ‖,橫條擁有相同‖qᵢ‖²,直條擁有相同‖qⱼ‖²。嘗試消去它們:每個橫條減去橫條總平均,然後每個直條減去直條總平均。不失一般性,假設原點位於數據中心,qᵢ總和是0,如此便剩下-2‖qᵢ ∙ qⱼ‖這一項。全部元素再除以-2,得到‖qᵢ ∙ qⱼ‖。最終形成了矩陣YᵀY,好算多了!簡言之:距離平方矩陣Mᵢⱼ²,行列中心化,除以-2,變成矩陣YᵀY。

問題變成:給定YᵀY,求得Y。

仿照PCA的手法,YᵀY實施特徵分解,得到答案Y = √ΛEᵀ。保留大的、捨棄小的特徵值和特徵向量,以符合新維度,同時也讓誤差最小。

YᵀY = EΛEᵀ = E√Λ√ΛEᵀ = (E√Λ)(√ΛEᵀ) = (√ΛEᵀ)ᵀ(√ΛEᵀ)

還有另一種手法,YᵀY實施Cholesky分解(對稱矩陣的LU分解),得到答案Y = Lᵀ。然而降維時誤差大,沒人這樣做。

YᵀY = LLᵀ

演算法(Random Projection)

https://en.wikipedia.org/wiki/Random_projection
https://en.wikipedia.org/wiki/Johnson–Lindenstrauss_lemma
https://homes.cs.washington.edu/~thickstn/docs/fast_jlt.pdf

線性函數的本質是旋轉、鏡射、伸縮座標軸,數據的相對位置不變。隨便找個線性函數做為嵌入函數,都能大致保留數據的相對位置,但是會大幅改變數據的相對距離。

而隨機導致分散,分散導致正交,正交導致距離不變。

絕聖棄智,民利百倍。

Manifold Embedding(Manifold Learning)

更換數據所在空間,從流形換成N維空間。

對付流形的手法主要是:一、取樣。二、建立鄰居關係(不見得是平面圖、三角剖分、網格)。三、鄰居關係盡量保持不變。

缺點是矩陣運算時間複雜度O(N³)太高。

Isomap: All Pair Shortest Path + Metric Multidimensional Scaling
Locally Linear Embedding: Weighted Average of Neighbor
Self-organizing Map: 1-layer Neural Network
Auto-encoder: Neural Network
https://en.wikipedia.org/wiki/Nonlinear_dimensionality_reduction
http://www.cs.utah.edu/~piyush/teaching/spectral_dimred.pdf
http://www.cad.zju.edu.cn/reports/流形学习.pdf
http://web.stanford.edu/class/ee378b/lect-9.pdf
http://www.cs.cmu.edu/~liuy/distlearn.htm

演算法(Stochastic Neighbor Embedding)

http://www.datakit.cn/blog/2017/02/05/t_sne_full.html

嵌入時,盡量保持內插函數相同。

數據添加一個新維度做為函數值,一律是1。原數據X、新數據Y,分別實施IDW Interpolation。令兩個內插函數的距離盡量小,距離採用KL Divergence。

距離倒數函數改成常態分布,即是此演算法,簡稱SNE。

距離倒數函數改成學生t分布,效果更好,簡稱t-SNE。

鐘形函數用於對付流形:超過一定距離,鄰居關係迅速減弱。

雖然細節有些出入,但是大致如此。

演算法(Locally Linear Embedding)

https://www.cnblogs.com/xbinworld/archive/2012/07/09/LLE.html

嵌入時,採用垂直投影?

變換矩陣,鄰居的加權平均,彷彿數據投影到座標軸。

多了鄰居概念的PCA?

演算法(Laplacian Eigenmap)

https://blog.csdn.net/jwh_bupt/article/details/8945083

嵌入時,盡量保持梯散相同。

鄰居矩陣A與梯散矩陣L:一、K個最近鄰居。二、heat kernel。

Laplace Smoothing改寫成Dirichlet Energy。

minimize xᵀLx   subject to ‖x‖² = 1

求出Laplace Matrix,找到第二小的特徵值暨特徵向量,特徵向量長度為1。(最小的特徵值0缺乏討論意義。)

原始論文將限制條件‖x‖² = xᵀx = 1改成了xᵀDx = 1。

備忘

從最佳化的觀點,xᵀLx的最小值是0,位於x=1。必須追加限制條件,以避開x=1。

原始論文追加兩個限制條件xᵀDx = 1與xᵀD1 = 0。xᵀDx = 1原理不明。xᵀD1 = 0讓加權平均值是零,xᵢ的權重是deg(xᵢ)。

以二維為例,xᵀDx = 1是橢圓,xᵀD1 = 0是過原點直線。單獨追加其中一個,答案仍是x=1。同時追加兩者,答案位於兩者交會之處。乍看像是隨手一指當作答案。

以下我另外提出一種更加鬼扯的解釋方式。

A可以正規化:每條邊的權重,除以端點的邊數開根號。如此一來,不會因為鄰邊太多,而導致梯散太大。

Asym = D-1/2 A D-1/2
Lsym = D-1/2 L D-1/2

找到一維數據,維度的共相關數總和越大越好(不是平方和),讓數據互相遠離。A用來設定哪些數據配對互相遠離。

maximize xᵀAx   subject to ‖x‖² = 1

預先將x的中心挪至原點。

maximize xᵀAx   subject to ‖x‖² = 1 and ∑xᵢ = 0

各點邊數不一致,誤差不平衡。除以度數,正規化。

maximize xᵀ√D⁻¹A√D⁻¹x   subject to xᵀx = 1 and xᵀ1 = 0
maximize yᵀAy           subject to yᵀDy = 1 and yᵀ√D1 = 0  (y := √D⁻¹x)

進一步改成Laplace Matrix。

minimize -yᵀAy          subject to yᵀDy = 1 and yᵀ√D1 = 0
minimize yᵀDy - yᵀAy    subject to yᵀDy = 1 and yᵀ√D1 = 0
minimize yᵀLy           subject to yᵀDy = 1 and yᵀ√D1 = 0

省略根號,如此一來就跟原始論文一樣。

minimize yᵀLy           subject to yᵀDy = 1 and yᵀD1 = 0

限制條件併入目標函數,化作廣義特徵值。

minimize yᵀLy - λ (yᵀDy - 1)   subject to yᵀD1 = 0
solve Ly - λDy = 0             subject to yᵀD1 = 0
solve Ly = λDy                 subject to yᵀD1 = 0

Deformation(Under Construction!)

Deformation

點邊面。

平面圖、球面圖,則可以討論面。

laplace operator 點是主角。

延伸應用

1. Graph Drawing: Stress Majorization。無向圖畫成圖片。
2. Fluid Simulation: Position-based Dynamics。調整水滴位置。
3. Structural Optimization: Topology Optimization。調整模型結構。
4. Molecular Modeling: Energy Minimization。調整化學結構。

無限制

min sum ‖(yi - yj) - (xi - xj)‖²

化作Laplace Matrix。各座標軸分開計算。

答案是所有點皆相等,縮成一點,缺乏討論意義。因此必須追加限制條件,實施約束最佳化。

精確位置限制

挑選一些點,設定明確座標。

解的一部分變成了定值,演算法只能使用Normal Equation系列、Jacobi Method系列,無法使用Conjugated Gradient Method系列。

相對位置限制、位置限制

min sum sum wij ‖xi - xj - vij‖² + sum ui ‖xi - pi‖²
min sum (Dx-v)ᵀW(Dx-v) + (x-p)ᵀU(x-p)
W and U are diagonal matrix
D is the discrete derivative operator
solve Lx = b such that L = DᵀWD + U and b = DᵀWv + Up

化作對稱對角優勢矩陣。

嚴格的說法是化作M-Matrix:Laplace Matrix的對角線元素各自加上任意正數或零。不過M-Matrix缺乏實務上的高速演算法。

相對距離限制

min sum kij (‖xi - xj‖ - dij)²
k是彈性係數 d是彈簧長度

似乎沒有公式解,只能使用數值模擬。請見本站文件「Solid Simulation」。

演算法是Euler Method:每個點分別求得加速度、更新速度、更新座標。

位置變化限制、速度變化限制、能量總和限制

min sum ‖xi' - xi‖² + ‖vi' - vi‖²  subject to c(x) = 0
x x'是先後位置,v v'是先後速度。c(x)是能量。

《FEPR: Fast Energy Projection for Real-Time Simulation of Deformable Objects》

二次函數有著高速演算法:Newton-Lagrange Method與Projected Gradient Method兩者合體。

乘子法-->微分方程組--->手動微分--->係數和乘數取出來變向量--->解大矩陣
qk+1 = argmin 1/2 ‖q - qₖ‖²D  subject to c(qₖ) + ∇c(qₖ)ᵀ(q - qₖ) = 0
         q
solve λk+1: [∇c(qₖ)ᵀ D⁻¹ ∇c(qₖ)] λk+1 = c(qₖ)
solve qk+1: qk+1 = qk - D⁻¹ ∇c(qk) λk+1

剛體變換限制、相似變換限制

min sum ‖(yi - yj) - T(xi - xj)‖²

Alignment + Deformation。

剛體變換是位移、旋轉。相似變換是位移、整體縮放、旋轉。

主角是點、主角是面,兩者式子不同。此處主角是點。

演算法是Block Coordinate Descent:固定R求得y,固定y求得R,兩者不斷輪流。換句話說,Procrustes Analysis與Laplace Equation不斷輪流。

Laplace Matrix是定值、是對稱矩陣,可以預先做Cholesky Decomposition,以便快速求解。

重心是定值,可以預先處理位移,將重心移到原點。

註:此演算法有人稱作Global-Local Method。

《As-Rigid-As-Possible Surface Modeling》,剛體,座標下降法,主角是點。

《Least Squares Conformal Maps for Automatic Texture Atlas Generation》,相似,Normal Equation,主角是面。

Manifold Deformation

請參考本站文件「Model」。

Smoothing(Under Construction!)

Smoothing

「平滑」。讓多筆數據的連線更加柔順。

Laplace Smoothing

梯散有兩種觀點。

加總差異:自身減鄰點,乘上權重,通通加總。
平均差異:自身減「鄰點的加權平均值」,乘以鄰點數。

每一點移動到「鄰點的平均值」。

每一點都正規化,預先除以鄰點數。避免點多而走遠。

x' = D⁻¹Ax = x - D⁻¹Lx   正規化
x' = x + laplace(x)      fixed point iteration
x' = x - Lx              因為圖論專家給他變號了

權重可自訂。經典的設定方式,請參考本站文件「Weighted Average Coordinates」。

Separation(Under Construction!)

Separation

「分離」。找到組成數據的關鍵因子。

例如麥克風錄到一段演奏,嘗試分離出每種樂器的聲音。

例如相機拍到一個場景,嘗試分離出光線的來源與強度。

演算法(Independent Component Analysis)

假定數據是關鍵因子的加權平均值。

演算法(Wavelet Analysis)

自訂特殊數據,做為關鍵因子。

Quantization(Under Construction!)

Quantization(Vector Quantization)

「量化」。找到一個函數,讓數據經過函數變換之後,保留數值的重點,刪除數值的細節。

簡易的量化是四捨五入、無條件捨去、無條件進入。進階的量化是區分數量級。經典的量化是以圖表相較並沒有明顯差異

如果預先知道數據大致分布,還可以自訂量化結果。亦可推廣到高維度。

量化只有一筆關鍵因子,分離則有許多筆關鍵因子。

演算法(Independent Component Analysis)

實施ICA分解版本,只做旋轉、略過縮放、略過位移。數據分解成加權平均值的格式,取權重最大者做為量化結果。

缺點是量化種類無法超過數據維度。

演算法(Cluster Analysis)

實施分群,以群集中心做為量化結果;再實施分類,以分界線決定量化結果。

https://charlesmartin14.wordpress.com/2012/10/09/spectral-clustering/
https://en.wikipedia.org/wiki/Spectral_clustering

演算法(Location-Allocation Analysis)

分群時,直線距離改成了其他指標。

Correlation(Under Construction!)

Correlation

兩堆數據的相關性。

演算法(Canonical Correlation Analysis)

https://www.ijcai.org/Proceedings/09/Papers/207.pdf

X投影到a(視作一物),Y投影到b(視作一物),最大化相關係數(兩物點積,而非兩兩交叉相乘)。

                                cov(aᵀX,bᵀY)                  
max corr(aᵀX,bᵀY) = max --------------------------------------
a,b                 a,b  sqrt(cov(aᵀX,aᵀX)) sqrt(cov(bᵀY,bᵀY))
                                 aᵀXYᵀb            
                  = max ---------------------------
                    a,b  sqrt(aᵀXXᵀa) sqrt(bᵀYYᵀb) 

a推廣成許多向量(一組基底),b推廣成許多個向量(一組基底)。

max tr(AᵀXYᵀB)   subject to AᵀXXᵀA = I, BᵀYYᵀB = I
A,B

演算法(Partial Least Squares Regression)

http://users.cecs.anu.edu.au/~kee/pls.pdf

維度的相關數矩陣XYᵀ投影到w,最大化變異數。

max tr(WᵀXYᵀYXᵀW)   subject to WᵀXXᵀW = I
 W

Power Iteration與Hotelling Deflation求得特徵值暨特徵向量。

Ranking(Under Construction!)

Ranking(Link Analysis)

節點排名。

演算法(PageRank)

演算法(HITS algorithm)

演算法(SpringRank)

min sum kij (‖xi - xj‖ - dij)²
k是彈性係數 d是彈簧長度
kij = dij = 1

Recommendation(Under Construction!)

Recommendation

點積是相關程度。收集兩兩相關程度,形成矩陣;分解矩陣,推敲出原始事物。

topic model = collaborative filtering
TF-IDF
word2vector
matrix completion
factorization machine
tensor decomposition http://people.csail.mit.edu/moitra/docs/Tensors.pdf
Laplacian Assignment