Create a Code Snippet by Logging In or Creating an Account
Cluster a correlation matrix

Wed, 04/08/2009 - 01:43 pm
This code demonstrates how you can cluster data using the correlation matrix.
// Create a correlation matrix from fake data, in this case with 2 embedded patterns and 

// relative noise of 0.25.  

FakeCorrelationMatrix(2,0.25)

// Cluster the correlation matrix with 2 expected patterns.  

// Picks starting cluster guesses randomly and for high values of noise it might give a poor 

// result.   The /SEED=(ticks) flag makes it give random initial guesses on each run.  

CorrelationClustering(2,Fake_Corr_Matrix)

// Display the data with cluster assignment, the original correlation matrix, and the clustered correlation matrix.  

Graph0(); Graph1(); Graph2()

// Create fake data and a correlation matrix from it.  

Function FakeCorrelationMatrix(num_patterns,relative_noise)

    Variable num_patterns // The number of patterns in the fake data.  

    Variable relative_noise // The standard deviation of the noise relative to the signal.  

    Variable length=50 // Each pattern should have this length.  

    Variable num_signals=30 // The number of cells/channels/signals in the fake data.  

    Make /o/n=(length,num_patterns) Base_Patterns=gnoise(1) // Make  a bunch of canonical patterns.  

    Make /o/n=(length,num_signals) Data // The matrix of fake data.  

    Make /o/n=(length,num_signals) Waterfall_Colors=q // Colors for the waterfall plot.  

    Make /o/n=(num_signals) Pattern_Values=floor(abs(enoise(num_patterns))) // The true pattern identities, chosen as an integer 0 to num_patterns-1.  

    Variable i

    for(i=0;i<num_signals;i+=1)

        Data[][i]=Base_Patterns[p][Pattern_Values[i]] // Each row will be a different pattern.   

        Data+=gnoise(relative_noise) // Add noise to everything.  

    endfor

    MatrixOp /O Fake_Corr_Matrix=syncCorrelation(Data) // Compute the covariance matrix.  

    MatrixOp /O Variances=varcols(Data)

    Fake_Corr_Matrix/=sqrt(Variances[p]*Variances[q]) // Convert to degree of correlation.  

    KillWaves /Z Base_Patterns,Variances // Cleanup.  

End

// Cluster a correlation matrix by swapping rows (and columns).  

Function CorrelationClustering(num_patterns,Corr_Matrix)

    Variable num_patterns // Number of patterns that you expect to find.  

    Wave Corr_Matrix // The correlation (or covariance) matrix.  

    KMeans /INIT=1 /NCLS=(num_patterns) /OUT=2 /SEED=(ticks) Corr_Matrix // K-Means clustering.  

    Duplicate /o Corr_Matrix Clustered_Matrix // Prepared the clustered correlation matrix.  

    Duplicate /o W_KMMembers Sorting_Index; Sorting_Index=p

    Sort W_KMMembers,Sorting_Index // Create a sorting index to use to swap out rows (and columns).  

    Clustered_Matrix=Corr_Matrix[Sorting_Index[p]][Sorting_Index[q]] // Shuffle rows (and columns).  

    Wave /Z Waterfall_Colors

    if(WaveExists(Waterfall_Colors))

        Wave W_KMMembers // The pattern number that each signal most represents (the K-Means clustering result).  

        Waterfall_Colors=W_KMMembers[q] // Color the waterfall plot according to the clustering result.  

    endif

    KillWaves /Z M_KMClasses,W_KMMembers // Cleanup.  

End

Window Graph0() : Graph

    PauseUpdate; Silent 1       // building window...

    Display /W=(509.25,50,924,358.25) as "Raw Correlation Matrix"

    AppendImage/T Fake_Corr_Matrix

    ModifyImage Fake_Corr_Matrix ctab= {-1,*,RedWhiteBlue,0}

    ModifyGraph margin(left)=14,margin(bottom)=14,margin(top)=14,margin(right)=14

    ModifyGraph mirror=2

    ModifyGraph nticks=4

    ModifyGraph minor=1

    ModifyGraph fSize=9

    ModifyGraph standoff=0

    ModifyGraph tkLblRot(left)=90

    ModifyGraph btLen=3

    ModifyGraph tlOffset=-2

    SetAxis/A/R left

EndMacro

Window Graph1() : Graph

    PauseUpdate; Silent 1       // building window...

    Display /W=(2.25,386.75,438.75,743.75) as "Clustered Correlation Matrix"

    AppendImage/T Clustered_Matrix

    ModifyImage Clustered_Matrix ctab= {-1,1,RedWhiteBlue,0}

    ModifyGraph margin(left)=14,margin(bottom)=14,margin(top)=14,margin(right)=14

    ModifyGraph mirror=2

    ModifyGraph nticks=3

    ModifyGraph minor=1

    ModifyGraph fSize=9

    ModifyGraph standoff=0

    ModifyGraph tkLblRot(left)=90

    ModifyGraph btLen=3

    ModifyGraph tlOffset=-2

    SetAxis/A/R left

EndMacro

Window Graph2() : Graph

    PauseUpdate; Silent 1       // building window...

    NewWaterfall /W=(-0.75,47,510.75,362.75)Data as "Data"

    ModifyWaterfall angle=45, axlen= 0.6, hidden= 0

    ModifyGraph negRGB=(0,0,65535)

    ModifyGraph zColor(Data)={Waterfall_Colors,*,*,Rainbow}

EndMacro