projects / epclust.git / blame
| Commit | Line | Data |
|---|---|---|
| ad642dc6 BA |
1 | ## File : 05_cluster2step.r |
| 2 | ## Description : | |
| 3 | ||
| 4 | rm(list = ls()) | |
| 5 | ||
| 6 | setwd("~/Documents/projects/2014_EDF-Orsay-Lyon2/codes/") | |
| 7 | source("~/Documents/projects/2014_EDF-Lyon2/R/01_StBr.r") # aux functs | |
| 8 | source('http://eric.univ-lyon2.fr/~jcugliari/codes/functional-clustering.r') | |
| 9 | library(cluster) | |
| 10 | ||
| 11 | ## 1. Read auxiliar data files #### | |
| 12 | ||
| 13 | identifiants <- read.table("identifs.txt")[ ,1] | |
| 14 | dates0 <- read.table("datesall.txt")[, 1] | |
| 15 | dates <- as.character(dates0[grep("2009", dates0)]) | |
| 16 | rm(dates0) | |
| 17 | ||
| 18 | n <- length(identifiants) | |
| 19 | p <- length(dates) | |
| 20 | ||
| 21 | synchros09 <- as.matrix(read.table("~/tmp/2009_synchros200.txt")) | |
| 22 | #synchros09 <- as.matrix(read.table("~/tmp/2009_synchros200RC.txt")) | |
| 23 | nas <- which(is.na(synchros09)[1, ]) # some 1/1/2009 are missing | |
| 24 | synchros09[nas, 1] <- colMeans(synchros09[2:4, nas]) | |
| 25 | ||
| 26 | #synchros09 <- t(synchros09) # series on rows, use if necessary | |
| 27 | ||
| 28 | ## 2. Compute contributions #### | |
| 29 | auxDWT <- t(apply(synchros09, 1, toDWT)) | |
| 30 | matcontrib <- t(apply(auxDWT, 1, contrib, rel = TRUE, logit = TRUE)) | |
| 31 | rm(auxDWT) | |
| 32 | ||
| 33 | ## 3. Transform data & compute CI #### | |
| 34 | ci <- CI(matcontrib[-201, ]) # last row has the total | |
| 35 | tdata <- ci$tdata; rownames(tdata) <- rownames(matcontrib)[-201] | |
| 36 | ||
| 37 | hc <- hclust(dist(tdata[, ci$selectv]), method = "ward.D") | |
| 38 | plot(hc) | |
| 39 | ||
| 40 | #clust <- cutree(hc, 2) | |
| 41 | ||
| 42 | for(K in 2:30){ | |
| 43 | #K <- 3 | |
| 44 | pamfit <- pam(tdata[-201, ci$selectv], k = K, stand = FALSE) | |
| 45 | ||
| 46 | #table(pamfit$clustering) | |
| 47 | ||
| 48 | SC <- matrix(0, ncol = p, nrow = K) | |
| 49 | ||
| 50 | clustfactor <- pamfit$clustering | |
| 51 | # for(k in 1:K){ | |
| 52 | # clustk <- which(clustfactor == k) | |
| 53 | # if(length(clustk) > 0) { | |
| 54 | # if(length(clustk) > 1) { | |
| 55 | # SCk <- colSums(synchros09[which(clustfactor == k), ]) | |
| 56 | # } else { | |
| 57 | # SCk <- synchros09[which(clustfactor == k), ] | |
| 58 | # } | |
| 59 | # SC[k, ] <- SC[k, ] + SCk | |
| 60 | # rm(SCk) | |
| 61 | # } | |
| 62 | #} | |
| 63 | ||
| 64 | write.table(clustfactor, file = paste0("~/tmp/clustfactorRC", K, ".txt")) | |
| 65 | #write.table(clustfactor, file = "~/tmp/clustfactor3.txt") | |
| 66 | } | |
| 67 | ||
| 68 | # Plots | |
| 69 | layout(1) | |
| 70 | matplot(t(SC)[48*10 + 1:(48*30), ], type = 'l', ylab = '',col = 1:3, lty = 1) | |
| 71 | matplot(t(SC)[48*100 + 1:(48*30), ], type = 'l', ylab = '', col = 1:3, lty = 1) | |
| 72 | ||
| 73 | ||
| 74 |