Functional Clustering

Introduction

• Functional clustering analysis
• Example codes
• Author: Taekyung Kim (PhD.), Kwangwwon University, Korea

Before you start

• Be aware that you should install gfortran for Mac OSX.
• If you have a new M1 Mac, go to https://github.com/fxcoudert/gfortran-for-macOS/releases to get a new gfortran package.
• You also need to install EMCluster package prior to use FClust.

Source

• The source codes in this example page mostly come from the documents in fdapace.

Sample Example

Load libraries

library(fdapace)
library(dplyr)
library(ggplot2)
library(gridExtra)
library(scales)

Load example data

data("medfly25")

Make Functional pricipal component analysis input

• IDs : identification
• tVec : time vector
• yVec : y vector

Flies <- MakeFPCAInputs(medfly25$ID,medfly25$Days,medfly25$nEggs)

Clustering

• See also, (FClust)
• cmethod : ‘EMCluster’ (default), ‘kCFC’
• default options:
  • methodMuCovEst = “smooth” : method to estimate [cross-sectional(default), smooth]
  • FVEthreshold = 0.90 : fraction-of-variance-explained by the components (0-1)
  • methodBwCov = ‘GCV’ : bandwidth choice method for the smoothed covariance function
  • methodBwMu = ‘GCV’ : bandwidth choice method for the mean function

Chen and Maitra (2015) Chiou and Li (2007)

newClust <- FClust(
  Flies$Ly,Flies$Lt,k=2,
  optnsFPCA = list(
    methodMuCovEst="smooth",
    userBwCov=2,
    FVEthreshold=0.90
  )
)

Observation

After observing data closely, we found that clusters are differentiated by the number of eggs that could be laid. The threshold value may be 25.

veryLowCount = ifelse(
  sapply(unique(medfly25$ID), function(u) 
    sum(medfly25$nEggs[medfly25$ID==u]))>=25,1,2)
veryLowCount

  [1] 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 2 1 1 1 1 1 2 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1
 [48] 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1
 [95] 1 1 1 1 1 1 2 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 2 1 1 1 1 2 1 1 1 1 2 1 1 1 2 2 1 1 1 1 1 1 1 1
[142] 1 2 1 1 2 2 1 1 1 1 1 1 2 1 1 1 1 1 1 1 1 1 1 1 1 1 1 2 1 1 2 1 1 1 1 1 1 1 1 1 1 1 2 2 2 1 1
[189] 1 1 1 2 1 1 1 1 1 1 1 1 1 2 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 2 1 1 1 1 1
[236] 1 1 1 1 1 2 1 1 1 2 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 2 1 1 2 1 1 1 1 1 1 1 1 1 1 1 1 1 1
[283] 1 1 1 1 1 1 1 1 1 1 1 1 2 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1
[330] 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 2 1 1 2 1 2 1 1 1 1 1 1 1 1 1 1
[377] 1 2 1 1 2 1 1 2 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 2 1 1 1 1 1 1 2 1 1 1 1 1 1 2 1 1 1 1 2
[424] 1 2 1 1 2 1 1 1 1 1 1 2 1 1 1 1 1 1 1 1 2 2 2 1 1 1 2 1 1 1 1 1 1 1 1 1 1 1 1 1 1 2 1 1 1 1 1
[471] 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 2 1 1 1 1 1 1 1 1 1 1 1 1 2
[518] 2 1 1 1 1 1 1 2 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 2 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 2 1 1 1 1
[565] 1 1 1 1 1 1 1 1 1 1 1 2 1 2 1 1 1 2 1 1 1 1 1 1 1 1 1 1 2 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1
[612] 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 2 1 1 1 1 1 1 1 2 1 1 1 1 2 1 1
[659] 1 1 1 2 1 1 1 2 1 2 1 1 1 1 2 1 1 1 1 1 1 1 1 1 2 1 1 1 1 1 1 1 1 1 1 1 1 1 1 2 1 1 1 2 1 1 1
[706] 1 1 1 1 1 1 1 1 1 1 1 1 2 1 2 1 1 2 1 1 1 1 1 1 1 2 1 2 1 1 2 1 1 1 1 1 2 1 1 1 1 1 1 1 1 1 1
[753] 1 2 1 1 1 1 1 1 1 1 1 1 1 2 1 2 1 1 1 1 1 1 2 2 1 1 1 1 1 1 1 1 2 1 1 1 2

Thus, the accuracy is calculated by:

sum((veryLowCount==newClust$cluster))/length(unique(medfly25$ID))

[1] 0.9961977

c1flies <- medfly25[newClust$cluster==1,]
c2flies <- medfly25[newClust$cluster==2,]

c1fliesObj <- MakeFPCAInputs(c1flies$ID,c1flies$Days,c1flies$nEggs)
c2fliesObj <- MakeFPCAInputs(c2flies$ID,c2flies$Days,c2flies$nEggs)

c1fliesFpca <- FPCA(c1fliesObj$Ly,c1fliesObj$Lt,list(methodMuCovEst="smooth"))
plot(c1fliesFpca)

c2fliesFpca <- FPCA(c2fliesObj$Ly,c2fliesObj$Lt,list(methodMuCovEst="smooth"))
plot(c2fliesFpca)

par(mfrow=c(1,2))
CreateModeOfVarPlot(c1fliesFpca,main="Cluster 1")
CreateModeOfVarPlot(c2fliesFpca,main="Cluster 2")

par(mfrow=c(1,1))