תרגול 12 - PCA and K-means

תקציר התיאוריה - PCA

PCA הוא אלגוריתם מאוד נפוץ אשר משמש במקומות רבים על מנת למצוא יצוג נוח יותר לוקטורים על סמך מדגם נתון. אחד השימושים העיקריים של האלגוריתם הינו בכדי לבצע הורדת מימד של הוקטורים. (לייצג וקטורים בעזרת וקטור ממימד נמוך יותר).

הגדרות

בעבור מדגם נתון D={x(i)}i=1N\mathcal{D}=\{\boldsymbol{x}^{(i)}\}_{i=1}^N של NN וקטורים באורך DD נגדיר את הגדלים הבאים:

  • הממוצע של המדגם: xˉ=1Ni=1Nx(i)\bar{\boldsymbol{x}}=\frac{1}{N}\sum_{i=1}^N \boldsymbol{x}^{(i)}.
  • מטריצת הדגימות:

    X=((x(1)μ)(x(2)μ)(x(N)μ))X=\begin{pmatrix} - & (\boldsymbol{x}^{(1)}-\boldsymbol{\mu})^{\top} & -\\ - & (\boldsymbol{x}^{(2)}-\boldsymbol{\mu})^{\top} & -\\ & \vdots & \\ - & (\boldsymbol{x}^{(N)}-\boldsymbol{\mu})^{\top} & -\\ \end{pmatrix}
  • הקווריאנס האמפירי של המדגם: P=XXP=X^{\top}X.

נתייחס לפירוק (ליכסון) הבא: P=UΛUP=U\Lambda U^{\top} כאשר UU היא מטריצה אורתונורמלית אשר העמודות שלה הם וקטורים עצמיים של PP:

U=(u1u2uD)U=\begin{pmatrix} | & | & & | \\ \boldsymbol{u}_1 & \boldsymbol{u}_2 & \dots & \boldsymbol{u}_D \\ | & | & & | \end{pmatrix}

ו Λ\Lambda היא מטריצה אלכסונית אשר מכילה את הערכים העצמיים של PP:

Λ=(λ1000λ2000λD)\Lambda=\begin{pmatrix} \lambda_1 & 0 & \dots & 0 \\ 0 & \lambda_2 & & 0 \\ \vdots & & \ddots & \vdots \\ 0 & 0 & \dots & \lambda_D \\ \end{pmatrix}

כך שהערך העצמי λj\lambda_j מתאים לוקטור העצמי uj\boldsymbol{u}_j והערכים העצמיים מסודרים מהגדול לקטן: λ1λ2λD\lambda_1\geq\lambda_2\geq\dots\geq\lambda_D.

הטרנספורמציה אותה מבצע PCA

PCA מייצר מתוך מדגם נתון D={x(i)}i=1N\mathcal{D}=\{\boldsymbol{x}^{(i)}\}_{i=1}^N טרנספורמציה אפינית (affine = linear + offset) אשר ממפה וקטור x\boldsymbol{x} באורך DD לוקטור z\boldsymbol{z} באורך KDK\leq D. כאשר KK הוא קבוע אשר נבחר מראש. הטרנפורמציה הינה:

z=T(xxˉ)\boldsymbol{z}=T^{\top}(\boldsymbol{x}-\bar{\boldsymbol{x}})

כאשר TT הינה מטריצה המכילה את KK העמודות הראשונות של UU (זאת אומרת הוקטורים העצמיים המתאימים ל KK הערכים העצמיים הגדולים ביותר).

האיברים של z\boldsymbol{z} נקראים הרכיבים הראשיים (principal components) של x\boldsymbol{x}.

פרשנות גיאומטרית

הפעולה שאותה מבצעת הטרנספורמציה הינה:

  1. להזיז את הנקודות של המדגם כך שהמרכז שלהם יהיה בראשית.
  2. הטלה של הנקודות המוזזות על תת-המרחב שמוגדרת על ידי הוקטורים {uj}\{\boldsymbol{u}_j\}.

מוטיבציה ראשונה: מקסימום שונות

תחת האילוץ ש TT הינה מטריצה בגודל D×KD\times K בעלת עמודות אורתו-נורמאליות, הבחירה הנוכחית של TT ממקסמת את הגודל:

1Ni=1Nz(i)22\frac{1}{N}\sum_{i=1}^N\lVert\boldsymbol{z}^{(i)}\rVert_2^2

אשר מכונה לרוב השונות של אוסף הוקטורים {z(i)}i=1N\{\boldsymbol{z}^{(i)}\}_{i=1}^N (בפועל זה שווה ל trace(ZZ)\text{trace}(Z^{\top}Z), כאשר ZZ מוגדרת בדומה ל XX רק עם הוקטורים z(i)\boldsymbol{z}^{(i)}).

מוטיבציה שניה: מזעור שגיאת השחזור הריבועית

נסתכל על זוג טרנספורמציות אפיניות כלליות מ x\boldsymbol{x} ל z\boldsymbol{z} באורך KK, ומ z\boldsymbol{z} ל x~\tilde{\boldsymbol{x}}:

z=Ax+bx~=Cz+d\begin{aligned} \boldsymbol{z}=A\boldsymbol{x}+\boldsymbol{b}\\ \tilde{\boldsymbol{x}}=C\boldsymbol{z}+\boldsymbol{d} \end{aligned}

נסמן את שגיאת השיחזור הריבועית באופן הבא: i=1N(x~(i)x(i))2\sum_{i=1}^N(\tilde{\boldsymbol{x}}^{(i)}-\boldsymbol{x}^{(i)})^2.

הטרנצפורמציות שימזערו את שגיאת השיחזור הריבועית הינם:

z=U(xμ)x~=UTz+μ\begin{aligned} \boldsymbol{z}=U(\boldsymbol{x}-\boldsymbol{\mu})\\ \tilde{\boldsymbol{x}}=U^T\boldsymbol{z}+\boldsymbol{\mu} \end{aligned}

תקציר התיאוריה - K-Means

K-Means הוא אלגוריתם אשכול אשר מנסה לחלק את הדגימות במדגם ל KK קבוצות על סמך המרחק בין הדגימות.

סימונים

  • KK - מספר האשכולות (גודל אשר נקבע מראש).
  • Ik\mathcal{I}_k - אוסף האינדקסים של האשכול ה-kk. לדוגמא: I5={3,6,9,13}\mathcal{I}_5=\left\lbrace3, 6, 9, 13\right\rbrace
  • Ik|\mathcal{I}_k| - גודל האשכול ה-kk (מספר הפרטים בקבוצה)
  • {Ik}k=1K\{\mathcal{I}_k\}_{k=1}^K - חלוקה מסויימת לאשכולות

בעיית האופטימיזציה

בהינתן מדגם D={x(i)}i=1N\mathcal{D}=\{\boldsymbol{x}^{(i)}\}_{i=1}^N, K-Means מנסה למצוא את החלוקה לאשכולות אשר תמזער את המרחק הריבועי הממוצע בין כל דגימה לכל שאר הדגימות שאיתו באותו האשכול. זאת אומרת, K-means מנסה לפתור את בעיית האופטימיזציה הבאה:

argmin{Ij}k=1K1Nk=1K12Iki,jIkx(j)x(i)22\underset{\{\mathcal{I}_j\}_{k=1}^K}{\arg\min}\frac{1}{N}\sum_{k=1}^K\frac{1}{2|\mathcal{I}_k|}\sum_{i,j\in\mathcal{I}_k}\lVert\boldsymbol{x}^{(j)}-\boldsymbol{x}^{(i)}\rVert_2^2

הבעיה השקולה

נגדיר את מרכז המסה של כל אשכול כממוצע של כל הוקטורים באשכול:

μk=1IkiIkx(i)\boldsymbol{\mu}_k=\frac{1}{|\mathcal{I}_k|}\sum_{i\in\mathcal{I}_k}\boldsymbol{x}^{(i)}

ניתן להראות כי בעיית האופטימיזציה המקורית, שקולה לבעיה של מיזעור המרחק הממוצע של הדגימות ממרכז המסה של האשכול:

argmin{Ij}k=1K1Nk=1KiIkx(i)μk22\underset{\{\mathcal{I}_j\}_{k=1}^K}{\arg\min}\frac{1}{N}\sum_{k=1}^K\sum_{i\in\mathcal{I}_k}\lVert\boldsymbol{x}^{(i)}-\boldsymbol{\mu}_k\rVert_2^2

האלגוריתם

K-mean הוא אלגוריתם חמדן אשר בכל פעם משייך מחדש את הדגימות ומעדכן את המרכזים.

האלגוריתם מאותחל בצעד t=0t=0 על ידי בחירה אקראית של KK מרכזי מסה: {μk}k=1K\{\mu_k\}_{k=1}^K.

בכל צעד tt מבצעים את שתי הפעולות הבאות:

  1. עדכון מחדש של החלוקה לאשכולות {Ik}k=1K\{\mathcal{I}_k\}_{k=1}^K כך שכל דגימה משוייכת למרכז המסה הקרוב עליה ביותר. כלומר אנו נשייך כל דגימה x\boldsymbol{x} לפי:

    k=argmink[1,K]xμk22k=\underset{k\in[1,K]}{\arg\min} \lVert\boldsymbol{x}-\boldsymbol{\mu}_k\rVert_2^2

    (במקרה של שני מרכזים במרחק זהה נבחר בזה בעל האינדקס הנמוך יותר).

  2. עדכון של מרכזי המסה המסה על פי:

    μk=1IkiIkx(i)\boldsymbol{\mu}_k=\frac{1}{|\mathcal{I}_k|}\sum_{i\in\mathcal{I}_k}\boldsymbol{x}^{(i)}

    (אם Ik=0|\mathcal{I}_k|=0 אז משאירים אותו ללא שינוי)

תנאי העצירה של האלגוריתם הינו כשהאשכולות מפסיקות להשתנות.

אחת הדרכים הנפוצות לאיתחול של {μk}k=1K\{\mu_k\}_{k=1}^K היא לבחור kk נקודות מתוך המדגם.

תכונות

  • מובטח כי פונקציית המטרה (סכום המרחקים מהממוצעים) תקטן בכל צעד.
  • מובטח כי האלגוריתם יעצר לאחר מספר סופי של צעדים.
  • לא מובטח כי האלגוריתם יתכנס לפתרון האופטימאלי, אם כי בפועל במרבית המקרים האלגוריתם מתכנס לפתרון אשר קרוב מאד לאופטימאלי.
  • אתחולים שונים יכולים להוביל לתוצאות שונות.

תרגיל 12.1 - PCA

עבוד מדגם נתון של וקטורים ב R2\mathbb{R}^2 חושבו וקטור הממוצע ומטריצת הקוואריאנס הבאים:

xˉ=(00)\bar{\boldsymbol{x}}=\begin{pmatrix}0\\0\end{pmatrix} P=(3226)P=\begin{pmatrix} 3 & 2 \\ 2 & 6 \end{pmatrix}

1) איזה מהוקטורים הבאים מייצג את הכיוון הראשון u1\boldsymbol{u}_1 במטריצת ההטלה של PCA?

15(21),12(11),15(12),\frac{1}{\sqrt{5}}\begin{pmatrix} -2 \\ 1 \end{pmatrix},\qquad \frac{1}{\sqrt{2}}\begin{pmatrix} 1 \\ 1 \end{pmatrix},\qquad \frac{1}{\sqrt{5}}\begin{pmatrix} 1 \\ 2 \end{pmatrix},\qquad

2) חשבו את שני ה principal componnents של x=(1,0)x=(1,0)^{\top}.

פתרון 12.1

1)

נשתמש בעובדה ש u1\boldsymbol{u}_1 צריך להיות וקטור עצמי של PP ולכן מקיים Pu1=λ1u1P\boldsymbol{u}_1=\lambda_1\boldsymbol{u}_1. נבדוק מי מהוקטורים הבאים מקיים זאת:

Pu1=15(3226)(21)=15(42)=2u1P\boldsymbol{u}_1 =\frac{1}{\sqrt{5}}\begin{pmatrix} 3 & 2 \\ 2 & 6 \end{pmatrix} \begin{pmatrix} -2 \\ 1 \end{pmatrix} =\frac{1}{\sqrt{5}}\begin{pmatrix} -4 \\ 2 \end{pmatrix} =2\boldsymbol{u}_1 Pu1=12(3226)(11)=12(58)αu1P\boldsymbol{u}_1 =\frac{1}{\sqrt{2}}\begin{pmatrix} 3 & 2 \\ 2 & 6 \end{pmatrix} \begin{pmatrix} 1 \\ 1 \end{pmatrix} =\frac{1}{\sqrt{2}}\begin{pmatrix} 5 \\ 8 \end{pmatrix} \neq \alpha\boldsymbol{u}_1 Pu1=15(3226)(12)=15(714)=7u1P\boldsymbol{u}_1 =\frac{1}{\sqrt{5}}\begin{pmatrix} 3 & 2 \\ 2 & 6 \end{pmatrix} \begin{pmatrix} 1 \\ 2 \end{pmatrix} =\frac{1}{\sqrt{5}}\begin{pmatrix} 7 \\ 14 \end{pmatrix} =7\boldsymbol{u}_1

מכאן שגם הוקטור הראשון וגם השלישי הם וקטורים עצמיים. הוקטור הראשון בהטלה של PCA יהיה השלישי שכן הוא מתאים לערך עצמי גדול יותר:

u1=15(12),\boldsymbol{u}_1= \frac{1}{\sqrt{5}}\begin{pmatrix} 1 \\ 2 \end{pmatrix},\qquad

2)

הרכיב העיקרי (principal componant) הראשון יהיה נתון על ידי:

z1=u1(xμ)=15(12)(10)=15z_1=\boldsymbol{u}_1^{\top}(\boldsymbol{x}-\mu) =\frac{1}{\sqrt{5}}\begin{pmatrix} 1 & 2 \end{pmatrix} \begin{pmatrix} 1 \\ 0 \end{pmatrix} =\frac{1}{\sqrt{5}}

והרכיב השני יהיה:

z2=u2(xμ)=15(21)(10)=25z_2=\boldsymbol{u}_2^{\top}(\boldsymbol{x}-\mu) =\frac{1}{\sqrt{5}}\begin{pmatrix} -2 & 1 \end{pmatrix} \begin{pmatrix} 1 \\ 0 \end{pmatrix} =\frac{-2}{\sqrt{5}}

בעבור PCA עם k=2k=2 נקבל:

z=15(1,2)\boldsymbol{z}=\frac{1}{\sqrt{5}}(1,-2)^{\top}

תרגיל 12.2

נתונות (1+3α)n\left(1+3\alpha\right)n נקודות שונות:

  • nn נקודות בקואורדינאטות A=(6,6)A=\left(-6,6\right)
  • αn\alpha n נקודות בכל אחת מהקואורדינאטות B=(6,6),C=(8,6),D=(1,6)B=\left(6,6\right),C=\left(8,6\right),D=\left(1,-6\right)

(הנקודות יושבות אחת על השניה בכל קואורדינטה, ומצויירות כעיגולים רק לצורך השרטוט). רוצים לבצע אשכול של הנקודות ל3 אשכולות בעזרת K-Means.

1) מאתחלים את המרכזים על ידי בחירה אקראית של 3 מתוך ארבעת הנקודות A,B,C,D. לאילו חלוקות יתכנס האלגוריתם בעבור כל אחת מארבעת האתחולים האפשריים.

2) מהו האשכול האופטימאלי (הממזער את פונקציית המטרה)? רשמו את הפתרון כתלות בפרמטר α\alpha. (ניתן להניח כי בפתרון האופטימאלי כל הנקודות שנמצאות באותו המקום משוייכות לאותו האשכול)

3) האם קיים אתחול אשר בעבורו האלגוריתם לא יתכנס לפתרון האופטימאלי שמצאתם בסעיף הקודם? הדגימו.

פתרון 12.2

1)

נחשב את תוצאת האלגוריתם בעבור כל אחת מארבעת האתחולים:

מרכזים ב A,B ו C:

  • שיוך התחלתי (0a): נקודות בA,B ו C ישוייכו למרכז אשר הנמצא עליהם, והנקודות בD ישוייכו למרכז שבB.
  • עדכון מרכזים (0b): המרכז שב B יזוז לאמצע הדרך שבין הנקודות B ו D.
  • עדכון אשכולות (1a): הנקודת שבB ישוייכו כעת למרכז שבC.
  • עדכון מרכזים (1b): המרכז שבין B ל D יזוז לD, והמרכז שבC יזוז למחצית הדרך שבין B לC.

מרכזים ב A,B ו D:

  • שיוך התחלתי (0a): נקודות בA,B ו D ישוייכו למרכז אשר נמצא עליהם, והנקודות בC ישוייכו למרכז שבB.
  • עדכון מרכזים (0b): המרכז שב B יזוז לאמצע הדרך שבין הנקודות B ו C.

מרכזים ב A,C ו D:

  • שיוך התחלתי (0a): נקודות בA,C ו D ישוייכו למרכז אשר נמצא עליהם, והנקודות בB ישוייכו למרכז שבC.
  • עדכון מרכזים (0b): המרכז שב C יזוז לאמצע הדרך שבין הנקודות B ו C.

מרכזים ב B,C ו D:

  • שיוך התחלתי (0a): נקודות בB,C ו D ישוייכו למרכז אשר נמצא עליהם, והנקודות בA ישוייכו למרכז שבB.
  • עדכון מרכזים (0b): המרכז שב B יזוז לנקודה שהיא המרכז של הנקודות A ו B. (משום שכמות הנקודות בשתי הקבוצות שונה, נקודה זו היא לא אמצע הדרך בניהם).

השלב הבא של עידכון האשכולות תלוי במיקום של המרכז החדש.

מקרה 1: הנקודות ב-B קרובות יותר למרכז החדש מאשר למרכז שב-C ולכן האלגוריתם מסתיים.

מקרה 2, המרכז החדש רחוק יותר לנקודה B מאשר הנקודה C, אזי הנקודות בB יהיו מושייכות כעת למרכז בנקודה C, והמשך האלגוריתם יהיה:

נמצא את התנאי על α\alpha שבעבורו מתרחש מקרה 2. נסמן ב μ1\boldsymbol{\mu}_1 את המרכז שבין A לB לאחר עדכון המרכזים הראשון. המיקום של μ1\boldsymbol{\mu}_1 נתון על ידי הממוצע המשוכלל של הקואורדיאנטות A ו B:

μ1=nA+αnB(1+α)n=(6x^1+6x^2)+α(6x^1+6x^2)1+α=α1α+16x^1+6x^2\boldsymbol{\mu}_1=\frac{n\vec{A}+\alpha n\vec{B}}{\left(1+\alpha\right)n}=\frac{\left(-6\hat{x}_1 + 6\hat{x}_2\right)+\alpha\left(6\hat{x}_1 + 6\hat{x}_2\right)}{1+\alpha}=\frac{\alpha-1}{\alpha+1}6\hat{x}_1 + 6\hat{x}_2

על מנת שיתרחש עידכון, על המרחק בין המרכז החדש לנקודה B צריך להיות גדול מ2:

(6x^1+6x^2)(α1α+16x^1+6x^2)>26α1α+16>2α1α+16<4α<5\begin{aligned} \left\lVert\left(6\hat{x}_1 + 6\hat{x}_2\right)-\left(\frac{\alpha-1}{\alpha+1}6\hat{x}_1 + 6\hat{x}_2\right)\right\rVert>2 \\ \Leftrightarrow 6-\frac{\alpha-1}{\alpha+1}6>2 \\ \Leftrightarrow \frac{\alpha-1}{\alpha+1}6<4 \\ \Leftrightarrow\alpha<5 \end{aligned}

2)

ב) אנו מועניינים למצוא את האשכול אשר מביא למינימום את הפונקציית המטרה הבאה:

argmin{Ij}k=1K1Nk=1KiIkx(i)μk22\underset{\{\mathcal{I}_j\}_{k=1}^K}{\arg\min}\frac{1}{N}\sum_{k=1}^K\sum_{i\in\mathcal{I}_k}\lVert\boldsymbol{x}^{(i)}-\boldsymbol{\mu}_k\rVert_2^2

נוכל לפסול פתרונות בהן ישנו אשכול ריק, משום שבמקרה זה נוכל לשייך אליו נקודות כלשהן על מנת להקטין את פונקציית המטרה. לכן הפתרון האופטימאלי חייב להיות אחד מששת האישכולים הבאים:

  • (A,B), (C), (D)
  • (A,C), (B), (D)
  • (A,D), (B), (C)
  • (B,C), (A), (D)
  • (B,D), (A), (C)
  • (C,D), (A), (B)

התרומה של האשכולות שמכילים קבוצה בודדת לפונקציית המטרה הינה 0, ולכן יש לחשב רק את התרומה של האשכול שמכיל שתי קבוצות של נקודות. למשל, עבור האשכול (A,B), (C), (D) נקבל:

argmin{Ij}k=1K1Nk=1KiIkx(i)μk22=n(66α1α+1)2+αn(66α1α+1)2=n36(α+1)2(4α2+4α)=144αnα+1\underset{\{\mathcal{I}_j\}_{k=1}^K}{\arg\min}\frac{1}{N}\sum_{k=1}^K\sum_{i\in\mathcal{I}_k}\lVert\boldsymbol{x}^{(i)}-\boldsymbol{\mu}_k\rVert_2^2 =n\left(-6-6\frac{\alpha-1}{\alpha+1}\right)^2 + \alpha n\left(6-6\frac{\alpha-1}{\alpha+1}\right)^2=n\cdot \frac{36}{\left(\alpha+1\right)^2}\left(4\alpha^2+4\alpha\right)=\frac{144\alpha n}{\alpha+1}

ועבור האשכול (B,C), (A), (D) נקבל:

argmin{Ij}k=1K1Nk=1KiIkx(i)μk22=αn(1)2+αn(1)2=2αn\underset{\{\mathcal{I}_j\}_{k=1}^K}{\arg\min}\frac{1}{N}\sum_{k=1}^K\sum_{i\in\mathcal{I}_k}\lVert\boldsymbol{x}^{(i)}-\boldsymbol{\mu}_k\rVert_2^2 =\alpha n\left(1\right)^2 + \alpha n\left(1\right)^2=2\alpha n

נחשב את הערך של פונקצייות המטרה בעבור כל אחד מששת האשכולים:

Clusters Objective
(A,B), (C), (D) 144αnα+1144\frac{\alpha n}{\alpha+1}
(A,C), (B), (D) 193αnα+1193\frac{\alpha n}{\alpha+1}
(A,D), (B), (C) 196αnα+1196\frac{\alpha n}{\alpha+1}
(B,C), (A), (D) 2αn2\alpha n
(B,D), (A), (C) 30.5αn30.5\alpha n
(C,D), (A), (B) 42.5αn42.5\alpha n

נשים לב כי הפתרון האופטימאלי יהיה חייב להיות (A,B),(C),(D) או (B,C),(A),(D) (משום שכל השאר בהכרח גדולים מהם). נבדוק בעבור אלו ערכים של α\alpha האשכול הראשון הינו האופטימאלי:

144αnα+1<2αnα>71\begin{aligned} 144\frac{\alpha n}{\alpha+1}<2\alpha n \\ \Leftrightarrow \alpha>71 \end{aligned}

אם כן, בעבור α>71\alpha>71 הפתרון האופטימאלי הינו (A,B),(C),(D) ובעבור α<71\alpha<71 הפתרון האופטימאלי הינו (B,C),(A),(D).

נסכם כי עבור אתחול המרכזים בנקודות B,C ו-D נקבל:

  • עבור α<5\alpha<5 האלגוריתם ישדך את B ו-C וזהו הפתרון האופטימאלי גלובלית.
  • עבור α>71\alpha>71 האלגוריתם ישדך את A ו-B וזה הפתרון האופטימאלי גלובלית.
  • עבור 5<α<715<\alpha<71 האלגוריתם ישדך את A ו-B אולם זהו אינו הפתרון הגלובלי.

נבדוק בעבור האתחולים מהסעיף הקודם, מהם המקרים שבהם האלגוריתם אינו מתכנס לפתרון האופטימאלי:

  • בעבור α>71\alpha>71 הפתרון האופטימאלי הינו (A,B),(C),(D), אך עבור 3 מתוך 4 האיחולים שבדקנו האלגוריתם התכנס לפתרון של (B,C),(A),(D).
  • בעבור α<71\alpha<71 הפתרון האופטימאלי הינו (B,C),(A),(D), אך במקרה של α>5\alpha>5 ואתחול של מרכזים ב B,C ו D מתקבל הפתרון של (A,B),(C),(D).

ג) כל המקרים שצויינו בסעיף הקודם. בנוסף, ניתן לדוגמא לאתחל שניים מתוך שלושת המרכזים בנקודות מאד רחוקות, ואז כל הנקודות ישוייכו למרכז השלישי.

חלק מעשי - מיקום חניונים בניו יורק

תזכורת: מדגם נסיעות המונית ב New York

נחזור למדגם של נסיעות מונית בניו-יורק בו השתמשנו בתרגולים הראשונים לחיזוי זמן הנסיעה. נציג את 10 הדגימות הראשונות במדגם (סה"כ במדגם זה 100,000 נסיעות).

passenger count trip distance payment type fare amount tip amount pickup easting pickup northing dropoff easting dropoff northing duration day of week day of month time of day
0 2 2.76806 2 9.5 0 586.997 4512.98 588.155 4515.18 11.5167 3 13 12.8019
1 1 3.21868 2 10 0 587.152 4512.92 584.85 4512.63 12.6667 6 16 20.9614
2 1 2.57494 1 7 2.49 587.005 4513.36 585.434 4513.17 5.51667 0 31 20.4128
3 1 0.965604 1 7.5 1.65 586.649 4511.73 586.672 4512.55 9.88333 1 25 13.0314
4 1 2.46229 1 7.5 1.66 586.967 4511.89 585.262 4511.76 8.68333 2 5 7.70333
5 5 1.56106 1 7.5 2.2 585.926 4512.88 585.169 4511.54 9.43333 3 20 20.6672
6 1 2.57494 1 8 1 586.731 4515.08 588.71 4514.21 7.95 5 8 23.8419
7 1 0.80467 2 5 0 585.345 4509.71 585.844 4509.55 4.95 5 29 15.8314
8 1 3.6532 1 10 1.1 585.422 4509.48 583.671 4507.74 11.0667 5 8 2.09833
9 6 1.62543 1 5.5 1.36 587.875 4514.93 587.701 4513.71 4.21667 3 13 21.7831

הבעיה: מציאת חניונים

חברת מוניות רוצה לשכור KK מגרשי חניה ברחבי העיר NYC בהם יוכלו לחכות המוניות שלה בין הנסיעות.

לשם כך היא מעוניינת לבחור באופן אופטימאלי את המיקומים של מגרשי החניות האלו כך שהמרחק הממוצע מנקודת הורדת הנוסע למרגש החניה הקרוב יהיה מינימאלי.

שדות רלוונטיים

הפעם נתמקד בשני השדות הבאים מהמדגם:

  • dropoff_easting - הקואורדינאטה האורכית (מזרח-מערב) של סיום הנסיעה
  • dropoff_northing - הקואורדינאטה הרוחבית (צפון-דרום) של סיום הנסיעה

(למתעניינים: הקואורדינאטות נתונות בUTM-WGS84, היחידות הן בקירוב קילומטר).

ויזואליזציה של נקודות ההורדה

הגדרה פורמאלית של הבעיה

נשתמש בסימונים הבאים:

  • x\mathbf{x} הוקטור האקראי של מיקום סיום הנסיעה
  • NN: מספר הנסיעות במדגם.
  • x(i)\boldsymbol{x}^{(i)} הוקטור של מיקום סיום הנסיעה ה ii.
  • ck\boldsymbol{c}_k: המיקום של מגרש החניה ה kk.

המטרה: למצוא את מיקומי החניונים האופטימאליים אשר ממזערים את הגודל הבא

{ck}k=1K=argmin{ck}k=1KE[minkxck2]\{c_k\}_{k=1}^{K*}= \underset{\{c_k\}_{k=1}^{K}}{\arg\min}\quad \mathbb{E}\left[\min_{k}\lVert\mathbf{x}-\boldsymbol{c}_k\rVert_2\right]

מכיוון שהפילוג האמיתי של x\mathbf{x} לא ידוע ננסה למזער את התוחלת האמפירית:

{ck}k=1K=argmin{ck}k=1K1Niminkx(i)ck2\{c_k\}_{k=1}^{K*}= \underset{\{c_k\}_{k=1}^{K}}{\arg\min}\quad \frac{1}{N}\sum_{i}\min_{k}\lVert\boldsymbol{x}^{(i)}-\boldsymbol{c}_k\rVert_2

נרשום את הבעיה על ידי חלוקת המדגם לאשכולות. נגדיר את האשכול Ik\mathcal{I}_k כאוסף של כל הנסיעות שהחניון ה kk הוא הקרוב ביותר לנקודת הסיום שלהן. באופן זה נוכל לרשום את בעיית האופטימיזציה באופן הבא:

{ck}k=1K=argmin{ck}k=1K1Nk=1KiIkx(i)ck2\{c_k\}_{k=1}^{K*}= \underset{\{c_k\}_{k=1}^{K}}{\arg\min}\quad \frac{1}{N}\sum_{k=1}^K\sum_{i\in\mathcal{I}_k}\lVert\boldsymbol{x}^{(i)}-\boldsymbol{c}_k\rVert_2

פתרון באמצעות K-Means

נשים לב כי הבעיה שקיבלנו דומה מאד לבעיה אותה K-Means מנסה לפתור, עם הבדל משמעותי אחד. K-Means ממזער את המרחק הריבועי הממוצע בעוד שאנו מחפשים למזער את המרחק האוקלידי. ישנם אלגוריתמים מורכבים יותר אשר פותרים את הבעיה שלנו, אך לבינתיים נשאר עם K-Means.

נציין שזהו מצב נפוץ שבו איננו מסוגלים לפתור בעיה מסויימת באופן ישיר אז אנו פותרים בעיה דומה לה בתקווה לקבל תוצאות מספקות, אך לא בהכרח אופטמאליות.

נשתמש באלגוריתם K-means על מנת לבחור את המיקום של 10 מגרשי חניה.

המרחק נסיעה הממוצע המתקבל הינו:

1Nk=1KiIkx(i)ck2=700m\frac{1}{N}\sum_{k=1}^K\sum_{i\in\mathcal{I}_k}\lVert\boldsymbol{x}^{(i)}-\boldsymbol{c}_k\rVert_2=700m

חושב לציין שהפתרון הזה הוא לא בהכרח הפתרון האופטימאלי משתי סיבות:

  1. K-Mean לא מבטיח התכנסות למינימום הגלובלי. דרך אחת לשפר את תוצאות האלגוריתם הינה להריץ אותו מספר פעמים עם איתחולים שונים.
  2. כפי שציינו קודם K-Mean ממזערת את השגיאה הריבועית הממוצעת. ניתן אם כן לשפר קלות את התוצאות על ידי שמירה על האשכולות אך תיקון המרכז לנקודה אשר ממזערת את המרחק עצמו.

הערה הנקודה אשר ממזערת את המרחק האוקלידי (בלי הריבוע) בינה לבין כל הנקודות באשכול נקראת החציון הגיאומטרי The Geometric Median (wiki). ניתן למצוא נקודה זו על ידי שימוש באלגוריתם המוכונה Weiszfeld's algorithm.

מציאת מספר החניונים האופטימאלי

עד כה השתמשנו ב10 חניונים, נרצה כעת לבחור גם מספר זה בצורה מיטבית. באופן כללי ככל שנגדיל את מספר החניונים מרחק הנסיעה לחניונים יקטן, אך מנגד התחזוקה של כל חניון עולה כסף.

נניח כי:

  1. עלות האחזקה של חניון הינה 10k$ לחודש.
  2. בכל חודש יהיו בדיוק 100k נסיעות.
  3. עלות הנסיעה של מונית בדרך לחניון הינה 3$ לקילומטר.

נרשום תחת הנחות אלו את העלות החודשית של אחזקת החניונים והנסיעה אליהם:

10K+1003E[minkxck2]10\cdot K+100\cdot3\cdot\mathbb{E}\left[\min_{k}\lVert\boldsymbol{x}-\boldsymbol{c}_k\rVert_2\right]

והמקבילה האמפירית:

10K+10031Nk=1KiIkx(i)ck210\cdot K+100\cdot3\cdot\frac{1}{N}\sum_{k=1}^K\sum_{i\in\mathcal{I}_k}\lVert\boldsymbol{x}^{(i)}-\boldsymbol{c}_k\rVert_2

מספר החניונים כ Hyper parameter

כעת עלינו לבצע אופטימיזציה גם על מספר החניונים וגם המיקום שלהם. ראינו כיצד ניתן למצוא פתרון בעבור KK נתון, אך אין לנו דרך פשוטה להכליל את זה ל KK כלשהו. נוכל לעבור על כל ערכי KK הרלוונטים, לפתור את הבעיה עבורם ולבסוף לקחת את הפתרון הטוב ביותר. KK הוא למעשה hyper-parameter של הבעיה.

נריץ את אלגוריתם ה K-Means בעבור כל ערך של K[1,25]K\in[1,25], נשרטט את עלות הנסיעה, עלות אחזקת החניונים והעלות הכוללת:

נקבל כי:

  • מספר החניונים האופטימאלי הינו: 12.
  • מרחק הנסיעה הממוצע יהיה 630 מ'.
  • העלות הכוללת תהיה 308.12k$ לחודש.