Nothing
## ----setup, include=FALSE-----------------------------------------------------
knitr::opts_chunk$set(echo = TRUE,eval=TRUE,fig.align="center",fig.width = 7,fig.height = 6)
## -----------------------------------------------------------------------------
library(ClustGeo)
data(estuary)
dat <- estuary$dat
head(dat)
D.geo <- estuary$D.geo
map <- estuary$map
## ----out.width="55%"----------------------------------------------------------
# description of 5 municipalities in the map
head(map@data[,4:8])
# plot of the municipalities
library(sp)
?"SpatialPolygonsDataFrame-class"
sp::plot(map, border="grey") # plot method
sel <- map$NOM_COMM%in% c("BORDEAUX", "ARCACHON", "ROYAN") # label of 3 municipalities
text(sp::coordinates(map)[sel,],
labels = map$NOM_COMM[sel])
# we check that the municipalities in map are the same than those in X
identical(as.vector(map$INSEE_COM),rownames(dat))
## ----out.width="80%"----------------------------------------------------------
D0 <- dist(dat) # the socio-economic distances
tree <- hclustgeo(D0)
plot(tree,hang = -1, label = FALSE,
xlab = "", sub = "",
main = "Ward dendrogram with D0 only")
rect.hclust(tree ,k = 5, border = c(4,5,3,2,1))
legend("topright", legend = paste("cluster",1:5),
fill=1:5,bty= "n", border = "white")
## ----out.width="70%"----------------------------------------------------------
# cut the dendrogram to get the partition in 5 clusters
P5 <- cutree(tree,5)
city_label <- as.vector(map$"NOM_COMM")
names(P5) <- city_label
plot(map, border = "grey", col = P5,
main = "Partition P5 obtained with D0 only")
legend("topleft", legend = paste("cluster",1:5),
fill = 1:5, bty = "n", border = "white")
## -----------------------------------------------------------------------------
# list of the municipalities in cluster 5
city_label[which(P5==5)]
## -----------------------------------------------------------------------------
D1 <- as.dist(D.geo) # the geographic distances between the municipalities
## -----------------------------------------------------------------------------
range.alpha <- seq(0,1,0.1)
K <- 5
cr <- choicealpha(D0, D1, range.alpha,
K, graph = FALSE)
cr$Q # proportion of explained inertia
## ----out.width="60%"----------------------------------------------------------
?plot.choicealpha
plot(cr)
## ----echo=FALSE,eval=FALSE----------------------------------------------------
# #postscript("plotQ-1.eps",width = 950,height = 489)
# plot(cr)
# #dev.off()
# #postscript("plotQnorm-1.eps")
# plot(cr,norm=TRUE)
# #dev.off()
## -----------------------------------------------------------------------------
tree <- hclustgeo(D0,D1,alpha=0.2)
P5bis <- cutree(tree,5)
## ----out.width="70%"----------------------------------------------------------
sp::plot(map, border = "grey", col = P5bis,
main = "Partition P5bis obtained with alpha=0.2
and geographical distances")
legend("topleft", legend=paste("cluster",1:5),
fill=1:5, bty="n",border="white")
## ----fig.height=4,fig.width=4,fig.align="center",warning=FALSE,message=FALSE----
library(spdep)
?poly2nb
list.nb <- poly2nb(map, row.names = rownames(dat)) #list of neighbours of each city
?nb2mat
A <- nb2mat(list.nb,style="B")
diag(A) <- 1
colnames(A) <- rownames(A) <- city_label
A[1:5,1:5]
## -----------------------------------------------------------------------------
D1 <- as.dist(1-A)
## ----out.width="70%"----------------------------------------------------------
range.alpha <- seq(0,1,0.1)
K <- 5
cr <- choicealpha(D0, D1, range.alpha,
K, graph=FALSE)
plot(cr)
## ----out.width='70%'----------------------------------------------------------
cr$Qnorm # normalized proportion of explained inertia
plot(cr, norm = TRUE)
## ----out.width='70%'----------------------------------------------------------
tree <- hclustgeo(D0, D1, alpha =0.2)
P5ter <- cutree(tree,5)
sp::plot(map, border="grey", col=P5ter,
main=" Partition P5ter obtained with
alpha=0.2 and neighborhood dissimilarities")
legend("topleft", legend=1:5, fill=1:5, col=P5ter)
Any scripts or data that you put into this service are public.
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.