[ppam-mpi.git] / code / script_clustering_by_pam.r

# Original R script (not required by the C program)


##################################################################
## File: script_clustering_by_pam.r
##
## Description:   Using PAM to clustering conso using pam
## Last modified: jan 2012 by JC
##
##################################################################

## Function: Converts a matrix of curves to the discret wavelet domain
toDWT <- function(x, filter.number = 6, family = "DaubLeAsymm"){ 
  x2   <- spline(x, n = 2^ceiling( log(length(x), 2) ),
                 method = 'natural')$y
  Dx2 <- wd(x2, family = family, filter.number = filter.number)$D
  return(Dx2)
}

## Function: Computes the absolute contribution of the wavelet's scale 
##           to the total total of the curve.
contrib <- function(x) { 
  J   <- log( length(x)+1, 2)
  nrj <- numeric(J)
  t0  <- 1
  t1  <- 0
  for( j in 1:J ) {
    t1     <- t1 + 2^(J-j)
    nrj[j] <- sqrt( sum( x[t0:t1]^2 ) )
    t0     <- t1 + 1
  }
  return(nrj)
}

## 1. Load libraries & data ####
library(wavethresh)
library(cluster)


# powerload is a matrix that contains on each line one observation
#           of length delta

## 2. DWT  ##################
delta <- ncol(powerload)
n     <- nrow(powerload)
Xdwt  <- t(apply(powerload, 1, toDWT)) # DWT over the lines of powerload.
Xnrj  <- t(apply(Xdwt, 1, contrib))    # Absolute contribution to the energy.

## 3. Cluster  ##############
K <- 8 # Number of clusters
Xnrj_dist <- dist(Xnrj)
Xnrj_pam  <- pam(as.dist(Xnrj_dist), K)
Commit	Line	Data
	1	# Original R script (not required by the C program)
	2
	3
	4	##################################################################
	5	## File: script_clustering_by_pam.r
	6	##
	7	## Description: Using PAM to clustering conso using pam
	8	## Last modified: jan 2012 by JC
	9	##
	10	##################################################################
	11
	12	## Function: Converts a matrix of curves to the discret wavelet domain
	13	toDWT <- function(x, filter.number = 6, family = "DaubLeAsymm"){
	14	x2 <- spline(x, n = 2^ceiling( log(length(x), 2) ),
	15	method = 'natural')$y
	16	Dx2 <- wd(x2, family = family, filter.number = filter.number)$D
	17	return(Dx2)
	18	}
	19
	20	## Function: Computes the absolute contribution of the wavelet's scale
	21	## to the total total of the curve.
	22	contrib <- function(x) {
	23	J <- log( length(x)+1, 2)
	24	nrj <- numeric(J)
	25	t0 <- 1
	26	t1 <- 0
	27	for( j in 1:J ) {
	28	t1 <- t1 + 2^(J-j)
	29	nrj[j] <- sqrt( sum( x[t0:t1]^2 ) )
	30	t0 <- t1 + 1
	31	}
	32	return(nrj)
	33	}
	34
	35	## 1. Load libraries & data ####
	36	library(wavethresh)
	37	library(cluster)
	38
	39
	40	# powerload is a matrix that contains on each line one observation
	41	# of length delta
	42
	43	## 2. DWT ##################
	44	delta <- ncol(powerload)
	45	n <- nrow(powerload)
	46	Xdwt <- t(apply(powerload, 1, toDWT)) # DWT over the lines of powerload.
	47	Xnrj <- t(apply(Xdwt, 1, contrib)) # Absolute contribution to the energy.
	48
	49	## 3. Cluster ##############
	50	K <- 8 # Number of clusters
	51	Xnrj_dist <- dist(Xnrj)
	52	Xnrj_pam <- pam(as.dist(Xnrj_dist), K)