Matematika pro geometrickou morfometrii (3)

Matematika pro geometrickou morfometrii (3) Ján Dupej ([email protected]) Laboratoř 3D zobrazovacích a analytických metod Katedra antropologie a genetiky člověka Přírodovědecká fakulta UK v Praze

Opakování • Prokrustovská transformace (analýza, registrace) • Idea

• Algoritmus • • • • •

Rozklad na dva problémy (dvojice a množina) Těžiště v počátku Jednotková velikost Minimalizace vzálenosti od průměru (ne vzájemnou) Suboptimální

Matematika pro geometrickou morfometrii (3)

• Nalezení rigidní transformace která minimializuje vzájemnou vzdálenost eliminuje rozdílnou polohu a velikost

2

Opakování • Prokrustovská transformace

• Výhody • 2D, 3D, rychlá, stabilní (malá změna vstupních dat = malá změna výsledku) • Jednotková velikost

• Nevýhody • Rozloží chybu mezi landmarky tak aby celková byla nejmenší • Všechny landmarky mají stejnou váhu – Pinocchio efekt


• Postupná aproximace průměrného tvaru • Zarovnání dvojice (jedinec na aproximaci průměrného tvaru)

3

Opakování • Stupeň volnosti = minimální počet čísel popisující tvar • Může jich být více, pak tvary tvoří zakřivený podprostor v tomto prostoru (např. po Prokrustovské analýze) • Statistické metody nefungují na zakřiveném prostoru • Projekce do tečné roviny pro malé vzdálenosti • Metody odvozené pro Kendallův prostor

• Vzdálenost mezi tvary • Proskrustovská • Částečná Prokrustovská • Úplná Prokrustovská – „indukuje“ Kendallův prostor


• Teorie tvaru

4

Zobrazení rozdílů


• Skupina nebo dvojice • Pro jednoduché konfigurace moho stačit překryté exempláře • Scatter plot

5

Zobrazení rozdílů


• Velikost, barva, orientace symbolů • Složitější ve 3D

6

Thin Plate Spline

• Znám posun landmarků • Co je mezi, spočítám interpolací

• 𝑇𝑃𝑆 𝒙 = 𝒂0 + 𝒂1 𝑥 + 𝒂2 𝑦 +

𝑛 𝑖=1 𝒘𝑖 𝑈

𝒑𝑖 − 𝒙

• 𝒙 = (𝑥, 𝑦) • Radial basis function: 𝑈 𝑡 = 𝑡 2 ln⁡(𝑡)

?


• Větší počet landmarků – mnoho šipek, nepřehledné • Zobrazení pomocí deformované mřížky • Vytvoření spojitého prostoru pomocí interpolace

7

TPS odvození Afinní část

Elastická část 𝑛

𝑤𝑖𝑥 𝑈 𝒑𝑖 − 𝒙 𝑖=1

Neznámé

• Znám: 𝑇𝑃𝑆𝑥 𝒑𝑗 = 𝒒𝑗𝑥 ⁡⁡𝑗 = 1, . . 𝑛 • 3 + 𝑛 neznámých, 𝑛 rovnic – chybí tři rovnice abychom měli jednoznačné řešení 𝑛

𝑛

𝑤𝑖𝑥 = 0 𝑖=1

𝑛

𝑤𝑖𝑥 𝑝𝑖𝑥 = 0 𝑖=1

𝑤𝑖𝑥 𝑝𝑖𝑦 = 0 𝑖=1


𝑇𝑃𝑆𝑥 𝒙 = 𝑎0𝑥 + 𝑎1𝑥 𝑥 + 𝑎2𝑥 𝑦 +

8

TPS odvození • Soustava rovnic 1, 𝑝𝑗𝑥 , 𝑝𝑗𝑦 , 𝑈1j , … , 𝑈𝑛𝑗 𝑛

𝑤𝑖𝑥 𝑈 ∥ 𝐩𝐢 − 𝐩𝐣 ∥ = 𝑞𝑗𝑥

𝑖=1

𝐰 = 𝑎0x , 𝑎1x , 𝑎2x , 𝑤1x , … , 𝑤𝑛𝑥

𝑛 daných podmínek 3 další podmínky

1 ⋮ 1 0 0 0

𝑝1x ⋮ 𝑝𝑛𝑥 0 0 0

𝑝1y ⋮ 𝑝𝑛𝑦 0 0 0

𝑈1j ⋮ 𝑈𝑛𝑗 1 𝑝1x 𝑝1y

… ⋮ … … … …

𝑈𝑛𝑗 𝑎0x 𝑞0x ⋮ 𝑎1x ⋮ 𝑈𝑛𝑛 𝑎2x 𝑞𝑛𝑥 = 0 1 𝑤1x ⋮ 0 𝑝𝑛𝑥 0 𝑝𝑛𝑦 𝑤𝑛𝑥


𝑇𝑃𝑆𝑥 𝐩𝑗 = 𝑎0x 1 + 𝑎1x 𝑝𝑗𝑥 + 𝑎2x 𝑝𝑗𝑦 +

9

TPS odvození • Soustava rovnic

• Získám potřebné parametry 𝒘 pro transformační funkci • Provedu vykreslení • Levá dolní část 𝑀−1 velikosti 𝑛 × 𝑛 se nazývá matice ohybové energie (Bending energy matrix)


𝑀𝒘 = 𝒙 𝑀−1 𝑀𝒘 = 𝑀−1 𝒙 𝒘 = 𝑀−1 𝒙

10


TPS ukázka • Morphome3cs

11


TPS ukázka

12

TPS - shrnutí • První použití visualizace D´Arcy Thompson, ruční kreslení

• Inspirace v matematice strojního inženýrství

• Interpolace definovaná všude • Může vést k falešným závěrům

• Důležité je husté pokrytí landmarky • Visualizace rozdílů jako deformací ve 2D • Lze snadno ve 3D, problém se zobrazováním • 2D řezy, 3D mřížky

• Landmarky se nesmějí „křížit“ – nestabilita


• Nekonečně tenký kovový plátek

13

Ukázka 1 • Jiné použití TPS – rektifikace obrazu • Lokální deformace má globální efekt Vznikající nestabilita

Vzdálená deformace Matematika pro geometrickou morfometrii (3)

Přitáhnout prsty k sobě

14


Ukázka 2

15

Warps • Parametrický model tvaru odvozený z TPS

• Vlastní vektory bending energy matrix 𝐸 • Setříděné nejvýraznější po dvou kolmé směry lokálních deformací

• Příklad 0 0 1 1 0.5 • 𝑃 = ⁡⁡ ⁡⁡ ⁡⁡ ⁡⁡ 0 1 1 0 0.5

𝑇

0 0 1 1 0.5 • 𝑄 = ⁡⁡ ⁡⁡ ⁡⁡ ⁡⁡ 0 1 1 0 0.75


• Principal warps, Partial warps

𝑇

16

Warps - příklad

0.4809 −0.2404 0.4809 −0.2404 −0.4809 −0.2404 0.4809 −0.2404 0.4809 −0.4809 𝐸 = 0.4809 −0.2404 0.4809 −0.2404 −0.4809 −0.2404 0.4809 −0.2404 0.4809 −0.4809 −0.4809 −0.4809 −0.4809 −0.4809 1.9236

• Principal warps tvoří komponenty deformační složky bez ohledu na cílové landmarky


• Na spočtení 𝑀−1 potřebuji pouze souřadnice vstupních landmarků, 𝑃

𝑣1 = −0.5000 0.5000 −0.5000 0.5000 0.0000 𝑣2 = 0.2236 0.2236 0.2236 0.2236 −0.8944

17

Warps - model • Konkrétní exemplář získáme z:

𝑄 = 𝑃𝐴 +

𝑣1 𝑘1𝑇 + 𝑣2 𝑘2𝑇 𝑖 𝑈 ∥ 𝐩𝐢 − 𝐱 ∥

• Koeficienty 𝑘 se spočítají ze znalosti vstupních 𝑃 a výstupních 𝑄 landmarků – z TPS koeficientů • Nazývají se partial warp scores a tvoří tvarové proměnné

• Je možné vizualizovat vliv jednotlivých warpů na tvar odděleně (nastavením ostatních 𝑘𝑖 = 0)


• Afinní transformace vstupních landmarků • Lineární kombinace warpů

18

Warps - dekompozice

+

⋅ 𝑘2x , 𝑘2y


⋅ 𝑘1x , 𝑘1y

19

Analýza konečných prvků

• Motivace v mechanice kontinua • Popis oblasti mezi landmarky jako celku – konečné uzavřené prvky většinou trojúhelníkového tvaru • Přirazení hodnoty popisující rozdíly v tvaru nebo velikosti

• FESA – Finite Element Scaling Analysis


• Spojitá deformace (TPS nebo jiná) nenabízí přímo žádnou míru (číslo) hodnotící deformaci nebo lokální změnu velikosti

20


Konečné prvky

21

Mechanika kontinua

• Strain tensor • 𝐹 = 𝜖𝑖,𝑗 =

1 2

𝛻𝒖 + 𝛻𝒖𝑇

• Rozklad na dvojici tenzorů velikosti a deformace • 𝐹 =𝑆+𝑇

• Redukce každého tensoru na jediný skalár (lokální rozdíl velikostí/tvarů)


• Je možné dodefinovat spojitou transformační funkci a analyzovat ji prostředky mechaniky kontinua (Cheverud and Richtmeier 1986)

22

FESA na trojúhelníkových sítích

• Tvar je dán poměrem velikostí vlastních čísel deformační matice 𝑥1 𝑃 = 𝑄𝐴 = 𝑥2 𝑥3

𝑥1 𝑦1 𝑦2 = 𝑥2 𝑥3 𝑦3

𝑦1 𝑦2 𝑎 𝑦3 𝑐

𝑏 𝑑


• Velikost je dána poměrem ploch odpovídajících si trojúhelníků, nebo kružnice a elipsy

23

Lineární model tvaru • Model – matematický popis jevu z reálného světa • Reálné příklady ho mohou více nebo méně splňovat

• Model tvaru lebky může být jeden „vhodný“ reprezentant • Splňuje ho pouze on sám, extrémy jsou daleko

• Může být průměrný reprezentant ze zkoumané populace • Nesplňuje ho nikdo, ale všichni jsou blízko

• Průměr + tolerance • Může popisovat i nemožné exempláře

• Parametrický model • Každý exemplář odpovídá číselné konfiguraci


• Příklad

24

Lineární model tvaru • Zohledňuje tendence vyskytující se v pozorované množině jedinců, vzorku

• Takový model dokáže popsat nekoněčně velkou populaci – generování jedinců • Matematický nástroj – analýza hlavních komponent • Principal component analysis (PCA) • Jeden exemplář se skladá z příspěvků různých komponent • Hledáme komponenty které přispívají nejvíc


• 𝑃1 , 𝑃2 , … 𝑃𝑛 → 𝐹 𝑘1 , 𝑘2 , … 𝑘𝑚

25

PCA • Data 𝒑𝑖 = 𝑝𝑖1 , 𝑝𝑖2 , … , 𝑝𝑖𝑛

• Průměr 𝒑0 • Hlavní komponenty (PC) 𝒌𝑖

• Komponenty jsou seřazeny dle míry přispění (zachycené variability) • Každého jedince lze zapsat jako součet průměru a lineární kombinace komponent • 𝒑𝑖 = 𝒑0 + • 𝑣𝑗 skóre

𝑛 𝑗=1 𝑣𝑗 𝒌𝑗


• Rozklad na 𝒑0 , 𝒌𝑖

26

𝒌2

𝒑0 Matematika pro geometrickou morfometrii (3)

PCA 𝒌1 𝒑𝑖

27

PCA • Jak najít množinu hlavních komponent?

• Kovariance popisuje „závislost“ dvou náhodných proměnných • Kovariační matice popisuje vzájemnou „závislost“ všech dvojic 1 𝑐𝑖𝑗 = 𝑁

𝑃𝑘𝑖 − 𝑃𝑖 𝑃𝑘𝑗 − 𝑃𝑗 𝑘

𝑐11 𝐶= ⋮ 𝑐1𝑛

⋯ 𝑐1𝑛 ⋱ ⋮ ⋯ 𝑐𝑛𝑛


• Komplikovaná teorie, vynecháme (tl, dr) • Odpovídá vlastním vektorům kovariační matice

28

PCA výpočet • Rozklad kovariační matice na vlastní vektory

• Můžu se rozhodnout kolik komponent chci zahrnou do svého modelu • Víc komponent – víc (potenciálně zbytečných) detailů, parametrů • Poměr součtu vlastních čísel vybraných komponent ku celkovému součtu odpovídá množství informace


• Dvojic vlastní vektor 𝒙, číslo 𝜆 je stejně jako souřadnic • Vlastní číslo určuje důležitost komponenty • 𝐶𝒙 = 𝜆𝒙

29

PCA – váha komponent

Procento informace


• Nejvíc informace je typicky obsaženo v prvních několika komponentách – vzorek není náhodný

30 komponenta

PCA – příklad 1 1.1 𝑃= 1.1 1.3

• Kovariační matice 2 1.2 1.2 1.1

1.5 4 2.5 3.8 1.5 3.5 1.2 3.7

0.016 −0.039 −0.026 −0.012 −0.039 0.18 −0.031 0.068 𝑐𝑜𝑣 𝑃 = −0.026 −0.031 0.32 0.022 −0.012 0.068 0.022 0.043

• Vlastní vektory 𝐯𝟏 = −0.8117

−0.3713 −0.1295

𝐯𝟐 = 0.5388 −0.2217

• Model 0.4318

λ1 = 0

0% informace

λ2 = 0.016

𝐯𝟑 = 0.2167 −0.8866 −0.1194 −0.3909

λ3 = 0.2114

2.87% informace 37.92% informace

𝐯𝟒 = 0.0618 0.1642

λ4 = 0.3301

59.21% informace

−0.0343 0.8120 −0.9838 −0.0377

31

• Průměr 𝐏𝟎 = 1.1250


• Data

1.3750 1.6750 3.7500

PCA – příklad • Zbývá dopočítat koeficienty (souřadnice, skóre) pro namodelování původních dat −0.2303 0.5125

𝑃0 = 1.1250

1.3750 1.6750 3.7500

𝑃0 + 𝐯𝟒 𝑘1,4 = 1.3463 1.7911 1.4746 𝑃0 + 𝐯𝟑 𝑘1,3 + 𝐯𝟒 𝑘1,4 = 1.3761 1.7990

1.5021 3.9572

𝑃0 + 𝐯𝟐 𝑘1,2 + 𝐯𝟑 𝑘1,3 + 𝐯𝟒 𝑘1,4 = 1.3375 1.7760 1.4099 𝑃0 + 𝐯𝟏 𝑘1,1 + 𝐯𝟐 𝑘1,2 + 𝐯𝟑 𝑘1,3 + 𝐯𝟒 𝑘1,4 = 1.0000 2.0000

3.7307

3.9743


𝐤 𝟏 = 0.4158 0.1040

1.5000 4.0000 = 𝑃1

• Zobrazení vybraných 2 až 3 koeficientů do grafu pro všechny exempláře (scatter plot)

32


PCA – demonstrace

33

PCA aplikace – Eigenfaces (1) • Yale faces dataset


• Registrované fotografie lidských tváří • 2D data převedeny na 1D „sploštěním“

34

Obrázek z http://www.cs.princeton.edu/~cdecoro/eigenfaces/

PCA aplikace – Eigenfaces (2) • PCA na vektorech tváří produkuje hlavní komponenty – „eigenfaces“


• Každou tvář je možné vyjádřit jako součet průměrné tváře a lineární kombinace „eigenfaců“

35


PCA aplikace – Eigenfaces (3) • Kolik hlavních komponent je potřeba sečíst?


• Ideálně všechny (někdy není možné) • Čím více, tím přesnější aproximace • Nejvýraznější příspěvky od komponent s nízkým indexem

36


PCA aplikace – Eigenfaces (4) • Do báze určené PC je možné převést data 𝒒⁡které nebyly zahrnuty do analýzy


• Takové data nemusí být přesně rekonstruovatelné • 𝑣𝑗 = 𝒒 − 𝒑0 ∙ 𝒌𝑗

37


PCA – závěr • PCA modeluje vztahy mezi landmarky pouze lineárně • Existují metody které dokáží zachytit nelineární vztahy

• Dobře separuje třídy – odlišnost se projevuje v prvních komponentách (shluková analýza) • Redukce dimenze – vytvoří stejný počet „nových“ proměnných, ale poslední nesou jen minimum informace, ty je možné zanedbat • Dopočítání chybějících dat


• Další užitečné použití

38

Matematika pro geometrickou morfometrii (3)

Recommend Documents