KMeans Algorithm

KMeans 알고리즘은 주어진 데이터들을 사전에 정해진 수만큼의 클러스터로 자동 분류하는 간단한 클러스터링 알고리즘의 일종이다. 이 알고리즘의 주된 아이디어는 각 클러스터의 중심을 어떻게 하면 잘 잡는가이다. 좋은 선택은 각각의 중심이 가능한 서로 멀리 떨어지게 잡아야 한다. 일단 클러스터 중심이 설정되면 각각의 데이터에 대해서 가장 가까이에 있는 클러스터에 할당시키면 된다. 다음으로 클러스터의 중심을 그 클러스터에 속한 데이터의 무게중심으로 재설정하고, 다시 각각의 데이터를 클러스터에 할당하는 작업을 반복한다. 이 과정은 클러스터의 중심이 더 이상 이동하지 않을 때까지 반복 시행하면 된다. 결과적으로 이 알고리즘은 클러스터의 중심에서 그 클러스터에 할당된 데이터와의 거리의 제곱을 그 클러스터의 모든 점과 전체 클러스터에 대해서 합한 값을 최소화시키는 클러스터 중심과 데이터의 분할을 찾는 것이 된다.

1. k개의 클러스터의 중심을 설정한다.
2. 설정된 중심을 가지고, 각 데이터를 분할한다.
3. 분할된 데이터를 이용하여서 클러스터의 중심을 재설정한다.
4. 중심이 변화가 없을 때까지 2,3 과정을 반복한다.

비용 함수의 관점에서 보면 KMeans 클러스터링은 다음 함수를

$$ \text{cost}=\sum_{k \in classes} \sum_{i \in class_k} \big| \text{data}_k(i) - \text{mean}_k \big|^2 $$

을 최소화시키는 분할을 하는 것이다. 그러나,  KMeans는 항상 최적의 비용 함숫값을 보장하지는 않는다. 그리고 알고리즘의 결과는 초기값의 설정에 의존하여서 다른 결과를 낼 수가 있다. 보다 큰 단점은 클러스터의 수를 사전에 정해주어야 한다는 사실인데, 잘못된 클러스터 개수의 설정은 성능에 크게 영향을 미친다. 자동의 클러스터의 수를 설정하면서 클러스터링을 할 수 있는 간단한 것으로는 isodata 알고리즘이 있다.

typedef double * vector_t ;
double kmeans_label(
             const vector_t *obs,             // observations;
             const int n_obs,                   /* in # of vectors */            
             const vector_t *mean,           // mean vectors;
             const int n_mean,                /* # of mean vectors */
             int* label,                           // cluster labeling;
             const int veclen)                 // dimension of feature;
{
    double sqerr=0;
    for (int i = 0; i < n_obs; i++) {
        vector_t c = obs[i];
        /* Get an estimate of best distance (min_dist) and (min_idx) */
        int min_idx = label[i];
        vector_t m = mean[min_idx];
        double min_sqdist = 0 ;
        for (int l = 0; l < veclen; l++) {
            double t = m[l] - c[l]; min_sqdist += t * t;
        }
        for (int j = 0; j < n_mean; j++) {
            vector_t m = mean[j];
            // save time by checking:: dist < min_sqdist=previous min_val;
            double sqdist = 0;
            for (l = 0; (l < veclen) && (sqdist < min_sqdist); l++) {
                double t = m[l] - c[l]; sqdist += t * t;
            }
            if (sqdist < min_sqdist) {
                min_sqdist = sqdist; min_idx = j;
            }
        }
        label[i] = min_idx;
        sqerr += min_sqdist;
    }
    return sqerr ;
};

int kmeans_update(
              const vector_t *obs, //observations;
              const int n_obs,       //# of observations;     
              vector_t *mean,       // mean vectors;
              const int n_mean,    // # of mean vectors;
              int *label,               // cluster labelling;
              const int veclen      // dimension of features;
              ) {
    std::vector<int> cnt(n_mean);
    for (int i = 0; i < n_mean; i++)
        for (int l = 0; l < veclen; l++) mean[i][l] = 0.0;
        
    for (i = 0; i < n_obs; i++) {
        ASSERT((0 <= label[i]) && (label[i] < n_mean));  
        vector_t m = mean[label[i]];
        cnt[label[i]]++;  
        vector_t c = obs[i];
        if (c == NULL)
            TRACE("No observations for %u, but expected up through %u", i, n_obs-1);
        for (int l = 0; l < veclen; l++) m[l] += c[l];
    }
    for (i = 0; i < n_mean; i++) {
        int j = cnt[i];
        if (j != 0) for (int l = 0; l < veclen; l++) mean[i][l] /= (double) j;
        else {
            TRACE("Empty cluster %u\n", i);
            return FALSE;
        }
    }
    return TRUE;
}

double kmeans(
       const vector_t *obs,     // observations;
       const int n_obs,          /* # of observations */
       vector_t *mean,           /* initial set of means */
       const int n_mean,        /* # of means (should be k_mean?) */
       int *label,                   /* cluster labeling*/
       const int veclen,         /* vector length of means and corpus */    
       const double min_conv_ratio,
       const int max_iter
    ) {
    double p_sqerr = DBL_MAX;
    int ret;
    double sqerr = kmeans_label( obs, n_obs,mean, n_mean, &label[0], veclen);
    double conv_ratio = (p_sqerr - sqerr) / p_sqerr;
    for (int i = 0; (i < max_iter) && (conv_ratio > min_conv_ratio); i++) {
        TRACE("kmeans iter [%u] %e ...\n", i, conv_ratio);
        ret = kmeans_update(obs, n_obs, mean, n_mean, &label[0], veclen );
        if (ret != TRUE)
            return (double)ret;
        p_sqerr = sqerr;
        sqerr = kmeans_label(  obs, n_obs,mean, n_mean, &label[0], veclen);
        conv_ratio = (p_sqerr - sqerr) / p_sqerr;
    }
    TRACE("kmeans n_iter %u sqerr %e conv_ratio %e\n", i, sqerr, conv_ratio);   
    return sqerr;
}