4. Domácí úkol: Anscombeho datasety

Úvod

V osmdesátých letech minulého století v článku [1] přišel britský matematik Francis John Anscombe se čtyřmi datasety tvořenými body v rovině (tzv. Anscombe’s quartet), které sice všechny mají shodné statistické ukazatele (jako průměry, střední odchylky, korelace, atp.), ale na první pohled vypadají poměrně radikálně odlišně. Z pohledu ukazatelů tedy "vypadají shodně", ale vizuálně jsou "velmi odlišné".

V samotném článku se Anscombe nezabývá způsobem, jakým tohoto výsledku dosáhl. Článek [2] z roku 2017 za zabývá přesně touto problematikou a ukazuje, že lze k datasetu bodů v rovině s danými statistickými ukazateli nalézt dataset více méně arbitrárního tvaru, který mý shodné uvedené statistické ukazatele (v rámci jisté přesnosti). Hlavním nástrojem je optimalizační technika simulovaného ochlazování.

V tomto posledním úkolu budeme implementovat Julia balíček řešící tento problém pomocí techniky volně inspirované článkem [2].

Odkazy na literaturu použité v tomto úvodu.

[1] Anscombe, F.J., Graphs in Statistical Analysis. The American Statistician 27, 1, 17–21. (1973)
[2] Matejka, J., Fitzmaurice, G., Same Stats, Different Graphs: Generating Datasets with Varied Appearance and Identical Statistics through Simulated Annealing. CHI '17: Proceedings of the 2017 CHI Conference on Human Factors in Computing Systems, 1290-1294. (2017)

Teoretický background

Mějme zadanou sadu bodů v rovině, kterou budeme reprezentovat dvourozměrným polem initial_ds (první index nechť rozlišuje $x$ -ovou a druhý $y$ -ovou souřadnici, druhý pak udává jednotlivé body, jde tedy o $2 \times n$ matici). Dále nechť je zadán cílový dataset target_ds, do jehož tvaru se budeme snažit počáteční dataset initial_ds deformovat při zachování jeho statistických ukazatelů.

Označme si jako current_ds aktuální stav datasetu (na počátku tedy initial_ds). Klíčovou snahou našeho algoritmu bude minimalizace funkce udávající součet minimálních vzdáleností jednotlivých bodů z current_ds od bodů z target_ds` takovým způsobem, aby statistické ukazatele nového datasetu `current_ds zůstaly shodné vzhledem k počátečnímu initial_ds, v zadané toleranci. Označme si tuto funkci symbolem fit (fitness).

Algorimus postupuje iterativně. Nyní popišme jeden krok algoritmu, kdy aktuální stav datasetu je current_ds a teplota (kladný parametr, který mezi kroky budeme měnit) je T.

Nejprve vytvořme přípustného kandidáta na nový stav datasetu, označme si ho test_ds. Stačí drobná úprava current_ds spočívající v drobné změně polohy jednoho náhodně vybraného bodu. Tato změna musí být ale taková, že zvolené statistické ukazatele datasetu test_ds jsou shodné s vlastnostmi datasetu initial_ds, v zadané přesnosti. (Tento požadavek lze případně relaxovat, ale požadavek na zachování vlastností by měl být alespoň vysoce pravděpodobný; například $95 %$ námi produkovaných změn splňuje tento požadavek).
Pokud nová hodnota fit(test_ds) je menší, než stará hodnota fit(current_ds), pak je změna přijata a test_ds nahradí current_ds. Nový dataset test_ds je ale přijat i v případě, kdy nedošlo k poklesu hodnoty funkce fit, a to s pravděpodobností exp(-(fit(test_ds) - fit(current_ds)) / T).

Tyto iterace opakujeme, dokud fit neklesne pod zadanou toleranci (hypotetické minimum je $0$ , tato situace ale nutně nemusí nastat), nebo dokud neprovedeme zadaný počet iterací max_iter.

Existují různé strategie pro měnění teploty, v článku [2] doporučují (a tím se můžeme řídit) nastavit teplotu v i-té iteraci dle vzorce max(minimal_T, initial_T * ( 1 - (i / max_iter)^2 )), kde initial_T je počáteční a minimal_T koncová zadaná teplota.

Případně můžete vyzkoušet další způsoby měnění teploty (při vysoké teplotě je pravděpodobnost změny nezmenšující fit větší, při nízké naopak nižší; snahou je aby iterace probíhaly jistý čas v obou režimech, nejprve prozkoumávání a poté zpřesňování).

Implementační pokyny

Balíček

Vytvořte Julia balíček/modul s názvem Anscombe, který bude možno lokálně instalovat příkazem

pkg> add "git@gitlab.fit.cvut.cz:BI-JUL/B241/USERNAME"#solution/04-anscombe

kde USERNAME je vaše uživatelské jméno. Pro přípravu struktury balíčku můžete použít na přednášce zmíněný balíček PkgTemplates.

Náš balíček bude exportovat metodu run!, která přijme počáteční dataset initial_ds, koncový dataset target_ds a dále parametry výpočtu jako keyword argumenty:

max_iter: počet iterací,
initial_T: počáteční teplota,
minimal_T: minimální teplota,
eps: tolerance odchylky mezi statistickými parametry (viz níže),
případně další dle vašeho výběru.

Vnitřní implementaci nechávám na vás, ale z popisu výše je přirozené implementovat dále minimálně metodu fit, výpočet kandidáta na nový dataset move_points a jeden krok algorimu perturb.

Hledání nového kandidáta

Doporučuji opakovaně vybírat jeden z bodů a lehce ho posunovat o vzdálenost R v náhodně vybraném směru, dokud maximum absolutní velikosti chyby mezi statistickými ukazately vstupního a aktuálního datasetu není menší než eps. Při opakování pokusu je vhodné hodnotu R změnit zpětným škalováním, tedy vynásobením R zvoleným číslem mezi $0$ a $1$ . Po přijetí změny pak naopak můžete R do další iterace opět zvětšit vynásobením zvoleným číslem větším než $1$ .

Statistické parametry

Pro jednoduchost použijme pouze čtyři parametry: průměrné hodnoty (mean) a standardní odchylky (std) $x$ -ových a $y$ -ových souřadnic bodů v datasetu. Změna je pak přijata, pokud absolutní hodnoty rozdílů těchto ukazatelů jsou menší než zadaná přesnost eps

Bonus

Pokuste se vytvořit animaci celého procesu (obrázek aktuálního stavu a jeho statistické parametry). Jako cílový stav vemte nějaký atraktivní tvar (například WolframAlpha dokáže poskytnout parametrizace zajímavých útvarů).

Pokuste se Váš balíček okrášlit o dokumentaci. Stačí, když ji půjde vygenerovat lokálně, v takovém případě mě na to upozorněte v MR.

Další poznámky a tipy

Při vlastním experimentu jsem datasety ve tvaru čtverce a kružnice o rozměrech cca $50 \times 50$ při natavení počáteční teploty na 0.4 a koncové na 0.01 a tolerancí změny statistických parametrů 1e-3 potřeboval k dosažení uspokojivého výsledku několik miliónů iterací. Vzdálenost při změně polohy jednoho bodu se pohybovala řádově kolem 0.01.

Testy

K balíčku přidejte alespoň několik základních testů vámi implementovaných funkcí/metod/struktur. Výslednou výpočetní funkci vzhledem k její netriviální výpočetní náročnosti explicitně netestuje. Případné testy (doporučuji minimálně otestovat výpočty vzdáleností a funkce fit) umístěte do adresáře test. Můžete se inspirovat předchozími úkoly, případně ukázkovými balíčky, které jsme procházeli ve výuce.

Při ladění vašeho programu je ideální jeho chování testovat na nějakém jednoduchém datasetu, například začít s pár desítkami bodů rozmístěných na hranici čtverce a jako cílový dataset vzít několik desítek bodů na kružnici (čtverec i kružnice se stejným těžištěm).

Řešení a odevzdání

Opět vytvořte větev odvozenou z větve assignment/04-anscombe a nezvěte ji solution/04-anscombe. Do solution/02-conway vložte své řešení editací souboru main.jl případně přidáním testů do složky test. Až budete se svým řešením spokojeni, vytvořte MR (to můžete i dříve, aspoň uvidíte výsledek testů, pokud je nespouštíte lokálně) a přiřaďte mě k němu jako assignee. Tímto aktem úkol odevzdáte.