Generalized Procrustes analysis

Standard

Previously we’ve seen how Procrustes analysis works: https://subokita.com/2014/04/07/procrustes-analysis/. Here’s a generalised version, that attempts to find a mean shape out from a set of points / shapes. 

For our test data, we randomly generated 6 set of points that are related to each other, each set of points is a version of first set of points that have already undergone translation, rotation, uniform scaling, and some sort of random jittering. Each of the points set are plotted in different color.

Screen Shot 2014-04-09 at 17.39.26

Basically we perform generalised Procrustes analysis, which is similar to original Procrustes analysis, but done repeatedly until convergence and a mean shape is found. The part of the image below that’s shaded shows the resulting mean shape:

Screen Shot 2014-04-09 at 17.39.33
Here’s the sample code of GPA, done in C++ using OpenCV:

/**
 * Recenter each matrix / set of points, by subtracting them with the mean / centroid
 */
vector Procrustes::recenter( const vector& X ) {
    vector result;
    for( Mat x : X ) {
        Scalar mu_x = cv::mean(x);
        Mat temp    = x - Mat(x.size(), x.type(), mu_x);
        result.push_back( temp );
    }
    return result;
}

/**
 * Normalize each set of points
 */
vector Procrustes::normalize( const vector& X ) {
    vector result;
    for( Mat x : X )
        result.push_back( x / cv::norm( x ) );
    return result;
}

/**
 * Perform rotation alignment between a set of points / shapes and a chosen mean shape
 * Points and mean shape should already be centered and normalize prior to the call of this function
 */
vector Procrustes::align( const vector& X, Mat& mean_shape ) {
    vector result;
    Mat w, u, vt;
    for( Mat x : X ){
        SVDecomp( mean_shape.reshape(1).t() * x.reshape(1) , w, u, vt);
        result.push_back( (x.reshape(1) * vt.t()) * u.t() );
    }
    return result;
}

/**
 * Perform a generalized Procrustes analysis to find the mean shape
 **/
vector Procrustes::generalizedProcrustes( std::vector& X, const int itol, const float ftol ) {
    /* Arbitrarily choose the first set of points as our mean shape */
    Mat mean_shape = X[0].reshape( 1 );
    
    int counter = 0;
    while( true ) {
        /* recenter, normalize, align */
        X = recenter( X );
        X = normalize( X );
        X = align( X, mean_shape );
        
        /* Find a new mean shape from all the set of points */
        Mat new_mean = Mat::zeros( mean_shape.size(), mean_shape.type() );
        for ( Mat x : X )
            new_mean += x;
        new_mean = new_mean / X.size();
        
        /* Perform the loop until convergence */
        float diff = norm( new_mean, mean_shape );
        if( diff <= ftol || counter > itol )
            break;
        
        mean_shape = new_mean;
        counter++;
    }
    
    /* Return the result as vector of points */
    vector result;
    mean_shape.reshape(2).copyTo( result );
    return result;
}

My version of generalized procrustes analysis code is adapted from https://github.com/catchagain/procrupy

You can view my example of generalised Procrustes analysis written using C++ and OpenCV over here: https://github.com/subokita/Sandbox/tree/master/Procrustes/Procrustes

Advertisements

Leave a Reply

Fill in your details below or click an icon to log in:

WordPress.com Logo

You are commenting using your WordPress.com account. Log Out / Change )

Twitter picture

You are commenting using your Twitter account. Log Out / Change )

Facebook photo

You are commenting using your Facebook account. Log Out / Change )

Google+ photo

You are commenting using your Google+ account. Log Out / Change )

Connecting to %s