Support : https://sigr2021.github.io/gwr/

library(tidyverse)
#> Error in library(tidyverse): there is no package called 'tidyverse'
library(sf)
#> Linking to GEOS 3.9.1, GDAL 3.3.0, PROJ 8.0.1
library(tmap)
library(plotly)
#> Error in library(plotly): there is no package called 'plotly'
library(gtsummary)
#> Error in library(gtsummary): there is no package called 'gtsummary'
library(GGally)
#> Error in library(GGally): there is no package called 'GGally'
library(GWmodel)
#> Error in library(GWmodel): there is no package called 'GWmodel'
library(spdep)
#> Error in library(spdep): there is no package called 'spdep'
library(reshape2)
#> Error in library(reshape2): there is no package called 'reshape2'

Méthode permettant d’analyser l’hétérogénéité spatiale

2 lois de géographie - dépendance spatiale : choses liées quand porche : auto-corrélation spatiale - hétérogénéité spatiale: les choses varient dans l’espace (moyenne ou variance, co-variance (les ))

Aujourd’hui étude de la covariance

Non-stationarité spatiale : instabilité dans l’espace des moyennes, variances et covariances.

Lien entre le prix de vente d’un bien immobilier et la surface généralement fort.

Exploration par l’analyse spatiale Méthode GWR : Geographical Weighted Regression : Régression géographiquement pondérée.

Cas d’études: étude du marché immobilier à Oléron.

4 étapes - chargement de la base - modèle de régression classique (moindes carrés) - GWR - extensions possibles - GWR multi-scalaire - classification sur les résidus de la GWR: régionalisation de l’espace

Import et préparation des données

Base DVF ouverte depuis 2019, réputée fiable

import de la base

data <- read_csv("https://files.data.gouv.fr/geo-dvf/latest/csv/2020/departements/17.csv.gz")
#> Error in read_csv("https://files.data.gouv.fr/geo-dvf/latest/csv/2020/departements/17.csv.gz"): could not find function "read_csv"
head(data) #Pour vérifier que l'import est correct
#>                                                                             
#> 1 function (..., list = character(), package = NULL, lib.loc = NULL,        
#> 2     verbose = getOption("verbose"), envir = .GlobalEnv, overwrite = TRUE) 
#> 3 {                                                                         
#> 4     fileExt <- function(x) {                                              
#> 5         db <- grepl("\\\\.[^.]+\\\\.(gz|bz2|xz)$", x)                     
#> 6         ans <- sub(".*\\\\.", "", x)

Filtrage

les communes
les ventes
uniquement les maisons
suppression des valeurs NA

dataOleron <- data %>% 
  filter(nom_commune %in% c("Dolus-d'Oléron",
                            "La Brée-les-Bains",
                            "Le Château-d'Oléron",
                            "Le Grand-Village-Plage",
                            "Saint-Denis-d'Oléron",
                            "Saint-Georges-d'Oléron",
                            "Saint-Pierre-d'Oléron",
                            "Saint-Trojan-les-Bains") & 
           nature_mutation == "Vente" & 
           type_local == "Maison" &
           !is.na(longitude) & 
           !is.na(surface_terrain) &
           !is.na(valeur_fonciere))
#> Error in attr(data, "tsp") <- c(start, end, frequency): object is not a matrix

Conversion en dataframe {sf}

dataSf <- dataOleron %>% 
  st_as_sf(coords = c("longitude","latitude"), 
           crs = 4326) # WGS84
#> Error in st_as_sf(., coords = c("longitude", "latitude"), crs = 4326): object 'dataOleron' not found

plot(st_geometry(dataSf))
#> Error in st_geometry(dataSf): object 'dataSf' not found

Import du fond de carte en shapefile

Données source: https://github.com/sigr2021/gwr

oleron <- st_read("../data-raw/oleron.shp")
#> Reading layer `oleron' from data source 
#>   `/Users/runner/work/notesSIGR2021/notesSIGR2021/data-raw/oleron.shp' 
#>   using driver `ESRI Shapefile'
#> Simple feature collection with 1 feature and 1 field
#> Geometry type: LINESTRING
#> Dimension:     XY
#> Bounding box:  xmin: -1.413073 ymin: 45.79904 xmax: -1.187825 ymax: 46.04829
#> Geodetic CRS:  WGS 84

tmap_mode("view")
#> tmap mode set to interactive viewing

## tmap mode set to interactive viewing

tm_shape(oleron) + 
  tm_lines(col = "black") + 
  tm_shape(dataSf) + 
  tm_dots(col = "red")
#> Error in as.list.environment(environment()): object 'dataSf' not found

Création d’une nouvelle variable : distance au littoral

Hypothèse : prix des maisons varient en fonction de la distance au littoral.

Distance au littoral : variable prédictive

dataSf$dist_litt <- st_distance(dataSf, oleron) %>% 
  as.numeric()
#> Error in st_crs(x): object 'dataSf' not found

Exploration des variables

Distribution de la variable dépendante (prix de vente)

plot_ly(dataSf, x = ~valeur_fonciere) %>% add_histogram()
#> Error in add_histogram(.): could not find function "add_histogram"

La forme de la distribution permet de choisir le modèle. C’est l’allure de la variable dépendante va aider au calibrage des données.

La distribution tire à droite, il faut utiliser une fonction log pour la rendre normale.

Suppression d’une valeur aberrante

dataSf <- dataSf %>% filter(valeur_fonciere > 1000)
#> Error in as.ts(x): object 'dataSf' not found

Distribution très dissymétrique

Transformation avec une fonction log() permet de normaliser la forme de la distribution

plot_ly(dataSf, x = ~log(valeur_fonciere)) %>% add_histogram()
#> Error in add_histogram(.): could not find function "add_histogram"

En régression, les valeurs abberrantes peuvent tirer les coefficients.

Distribution des variables indépendantes

a <- plot_ly(dataSf, x = ~log(dist_litt)) %>% add_histogram()
#> Error in add_histogram(.): could not find function "add_histogram"
b <- plot_ly(dataSf, x = ~log(surface_reelle_bati)) %>% add_histogram()
#> Error in add_histogram(.): could not find function "add_histogram"
c <- plot_ly(dataSf, x = ~log(surface_terrain)) %>% add_histogram()
#> Error in add_histogram(.): could not find function "add_histogram"
subplot(a,b,c)
#> Error in subplot(a, b, c): could not find function "subplot"

Fonction plot_ly permet d’avoir un graphique interactif.

En régression, les variables explicatives n’ont pas besoin d’être loggées mais c’est plus simple. Les modèles log-log, l’interprétation est plus simple (effet d’un pourcentage sur un autre pourcentage)

Poursuite de la préparation des données.

# Suppression des maisons vraisemblablement trop petites
dataSf <- dataSf %>% filter(surface_reelle_bati > 10)
#> Error in as.ts(x): object 'dataSf' not found

# Création des variables log (pour faciliter la carto par la suite)
dataSf$log_valeur_fonciere <- log(dataSf$valeur_fonciere)
#> Error in eval(expr, envir, enclos): object 'dataSf' not found
dataSf$log_dist_litt <- log(dataSf$dist_litt)
#> Error in eval(expr, envir, enclos): object 'dataSf' not found
dataSf$log_surface_reelle_bati <- log(dataSf$surface_reelle_bati)
#> Error in eval(expr, envir, enclos): object 'dataSf' not found
dataSf$log_surface_terrain <- log(dataSf$surface_terrain)
#> Error in eval(expr, envir, enclos): object 'dataSf' not found

Relations bivariées - formes fonctionelles

On regarde si les relations sont linéaires (pour pouvoir faire une régression linéaire)

Lien prix - distance au littoral


ggplot(dataSf, aes(x=log(dist_litt), y=log(valeur_fonciere))) + 
  geom_point() + geom_smooth()
#> Error in ggplot(dataSf, aes(x = log(dist_litt), y = log(valeur_fonciere))): could not find function "ggplot"

Relation assez plate mais pas non linéaire, légèrement négative (plus on s’éloigne du littoral plus prix baisse. Mais pas nécessairement, il y a de l’hétérogéinité spatiale).

Lien prix-surface réelle

ggplot(dataSf, aes(x=log(surface_reelle_bati), y=log(valeur_fonciere))) + 
  geom_point() + geom_smooth()
#> Error in ggplot(dataSf, aes(x = log(surface_reelle_bati), y = log(valeur_fonciere))): could not find function "ggplot"

Lien prix- surface terrain

ggplot(dataSf, aes(x=log(surface_terrain), y=log(valeur_fonciere))) + 
  geom_point() + geom_smooth()
#> Error in ggplot(dataSf, aes(x = log(surface_terrain), y = log(valeur_fonciere))): could not find function "ggplot"

Modèle log-log global (MCO)

Relation moyenne sur tout l’espace

Fonction lm : linear modele régression variable quantitative continue
Fonction glm : variable qualitative, modéle généralisé

mco <- lm(log_valeur_fonciere ~ log_dist_litt + log_surface_reelle_bati + log_surface_terrain, data = dataSf)
#> Error in is.data.frame(data): object 'dataSf' not found
mco %>%
  tbl_regression(intercept = TRUE) %>% add_vif() # problème (voir package gtsummary ?)
#> Error in add_vif(.): could not find function "add_vif"

VIF: Variance Inflation Factor

Si 0 n’est pas dans l’intervalle de confiance, cela veut dire que la relation est significative. Si 0 compris dans l’intervalle de confiance, la relation peut être nulle.

Interprétation si distance au littoral augmente de 1%, les prix diminuent de 0,08%.

summary(mco)
#> Error in summary(mco): object 'mco' not found

Représentation graphique des coefficients

ggcoef_model(mco)
#> Error in ggcoef_model(mco): could not find function "ggcoef_model"

surface bâti réel et surface terrain ont un effet positif significatif.

Cartographie des résidus

Résidus: différence entre la valeur y réelle et la valeur y du modèle.

On a souvent des résidus qui forment des pattern spatialisés

dataSf$resMco <- mco$residuals
#> Error in eval(expr, envir, enclos): object 'mco' not found

tm_shape(dataSf) + tm_dots(col = "resMco", style = "quantile")
#> Error in as.list.environment(environment()): object 'dataSf' not found

Couleur rouge : le modèle sous-estime le prix Couleur verte: le modèle sur-estime le prix

Modèle GWR

Principes généraux

MCO : un coefficient moyen par variable

GWR: coefficient moyen pour chaque maison avec le voisinage uniquement

Choix: distance du voisinage et pondération

Caractérisation du voisinage: Knn, semi-variogramme

GWR: cartographie des variations spatiales: endroits ou tel prédicteur joue et où il ne joue pas.

GWR, étape pour mieux comprendre le territoire et compléter son modèle.

Deux étapes

Définir le voisinage et ensuite le pondérer.

Caractériser le voisinage : distance fix, nbre de voisin, graphe…

Pondération: - binaire (voisin ou non) - pondération selon la distance (fonction inverse)

Il existe 2 packages R pour estimer la GWR.

GWModel est plus convivial mais non compatible avec SF

Première chose à faire : convertir l’objet sf en objet sp

dataSp <- as_Spatial(dataSf) # le package GWmodel n'est pas compatible avec 'sf'
#> Error in as_Spatial(dataSf): object 'dataSf' not found

Construction de la matrice de distances qui servira pour le calibrage. Donc cela peut devenir un peu long voir compliquer sur de gros échantillons.

matDist <- gw.dist(dp.locat = coordinates(dataSp))
#> Error in gw.dist(dp.locat = coordinates(dataSp)): could not find function "gw.dist"

Optimisation du voisinage

Paramètres important de GWR:

bw: fenêtre de voisinage
adaptative: distance (FALSE) ou nbre de voisin (TRUE)

Parfois problème de convergence si distance: manque de voisin, donc très grande distance. Nbre de voisin, on est certain d’avoir assez de voisin.

En général, on conseille d’utiliser un nombre de voisin si la répartition n’est pas régulière.

On cherche à minimiser l’AIC entre 2 fonctions pour choisir la meilleure

Fonction bi-carrée: passer un certain seuil la pondération est égale à 0

Comparaison de deux pondérations spatiales : exponentielle et bicarrée :

Calcul du nombre de voisin optimal pour chaque pondération spatiale

# Exponential
nNeigh.exp <- bw.gwr(data = dataSp, approach = "AICc", # fonction cible méthode de validation
                     kernel = "exponential", # fonction de noyau
                     adaptive = TRUE, # nbre de voisin
                     dMat = matDist,
                     formula = log_valeur_fonciere ~ log_dist_litt + log_surface_reelle_bati + log_surface_terrain)
#> Error in bw.gwr(data = dataSp, approach = "AICc", kernel = "exponential", : could not find function "bw.gwr"

Nombre de voisin optimal est égal à 25.

Valeur AICc: paramètre optimal pour le modèle

# Bisquare
nNeigh.bisq <- bw.gwr(data = dataSp, approach = "AICc", 
                      kernel = "bisquare", 
                      adaptive = TRUE, 
                      dMat = matDist, 
                      formula = log_valeur_fonciere ~ log_dist_litt + log_surface_reelle_bati + log_surface_terrain)
#> Error in bw.gwr(data = dataSp, approach = "AICc", kernel = "bisquare", : could not find function "bw.gwr"

Avec la fonction bisquare, Nombre de voisin optimal est égal à 151.

Estimation de la GWR

Avec pondération exponential

GWR.exp <- gwr.basic(data = dataSp, bw = nNeigh.exp, kernel = "exponential", adaptive = TRUE,  dMat = matDist, formula = log_valeur_fonciere ~ log_dist_litt + log_surface_reelle_bati + log_surface_terrain)
#> Error in gwr.basic(data = dataSp, bw = nNeigh.exp, kernel = "exponential", : could not find function "gwr.basic"

Avec pondération bisquare

GWR.bisq <- gwr.basic(data = dataSp, bw = nNeigh.bisq, kernel = "bisquare", 
                      adaptive = TRUE,  dMat = matDist, 
                      formula = log_valeur_fonciere ~ log_dist_litt + log_surface_reelle_bati + log_surface_terrain)
#> Error in gwr.basic(data = dataSp, bw = nNeigh.bisq, kernel = "bisquare", : could not find function "gwr.basic"

Comparaison des deux calibrations :

diagGwr <- cbind(
  rbind(nNeigh.exp,nNeigh.bisq),
  rbind(GWR.exp$GW.diagnostic$gw.R2,GWR.bisq$GW.diagnostic$gw.R2),
  rbind(GWR.exp$GW.diagnostic$AIC,GWR.bisq$GW.diagnostic$AIC)) %>% 
  `colnames<-`(c("Nb Voisins","R2","AIC")) %>% 
  `rownames<-`(c("EXPONENTIAL","BISQUARE"))
#> Error in rbind(nNeigh.exp, nNeigh.bisq): object 'nNeigh.exp' not found
diagGwr
#> Error in eval(expr, envir, enclos): object 'diagGwr' not found

Explication:

Exponentielle : maximise le R2 (explique 58% de la variation des prix) et minimise l’AIC : modèle plus performant

Bisquare: n’explique que 50% de la variation des prix, l’AIC est moins bon.

A partir de ces éléments, on va choisir la GWR exponentielle (plus performante).

Interprétation des résultats bruts de la GWR :

Interprétation des résultats de la GWR pas facile mais toujours intéressant.

GWR.exp
#> Error in eval(expr, envir, enclos): object 'GWR.exp' not found

Min et max distance : min négatif et max positif : inversion de signe par endroit.

Surface réelle: certains endroits la surface a un fort impact, d’autres moins

Cartographie des betas

pas de fonction toute prête pour cartographier les résidus

# Fonction de cartographie automatique des coefficients GWR et de la t-value
mapGWR <- function(spdf,var,var_TV,legend.title = "betas GWR",main.title, dot.size = 0.3) {
  tv <- spdf[abs(var_TV)>1.96,] # tv t-value: coeff régression / écart-type
  tm_shape(spdf) +
    tm_dots(var, title = legend.title, size = dot.size) +
    tm_shape(oleron) + tm_lines() +
    tm_shape(tv) + tm_dots(col="grey40") +
    tm_layout(title = main.title, legend.title.size =0.9, inner.margins = .15) 
}

t-value: coeff régression / écart-type

Certains endroits coefficients très élevés mais erreur standard aussi, la t-value le met en valeur.

Si t-value très basse : signal un coefficient

Planche cartographique des 3 variables

# 
tmap_mode("plot")
#> tmap mode set to plotting
a <- mapGWR(GWR.exp$SDF, var = "log_dist_litt",var_TV = GWR.exp$SDF$log_dist_litt_TV,
            main.title = "Distance au littoral")
#> Error in mapGWR(GWR.exp$SDF, var = "log_dist_litt", var_TV = GWR.exp$SDF$log_dist_litt_TV, : object 'GWR.exp' not found
b <- mapGWR(GWR.exp$SDF, var = "log_surface_reelle_bati",var_TV = GWR.exp$SDF$log_surface_reelle_bati_TV,
            main.title = "Surface bâtie")
#> Error in mapGWR(GWR.exp$SDF, var = "log_surface_reelle_bati", var_TV = GWR.exp$SDF$log_surface_reelle_bati_TV, : object 'GWR.exp' not found
c <- mapGWR(GWR.exp$SDF, var = "log_surface_terrain",var_TV = GWR.exp$SDF$log_surface_terrain_TV,
            main.title = "Surface terrain")
#> Error in mapGWR(GWR.exp$SDF, var = "log_surface_terrain", var_TV = GWR.exp$SDF$log_surface_terrain_TV, : object 'GWR.exp' not found

tmap_arrange(a,b,c)
#> Error in tmap_arrange(a, b, c): object 'a' not found

Variation de couleur: variation de la variable Si significatif : point au milieu du cercle.

Interprétation des cartes : - Surface bâtie : effet de la surface batie plus importante à chateau d’Oléron et à Saint Pierre d’Oléron (valeurs rouge) qu’ailleurs. Pour une étude très locale, cela ne vaudra pas le coup d’étudier cette variable - Surface terrain: Saint Denis - Distance au littoral: Saint Trojan, les prix vont beaucoup augmenter quand on se rapproche du littoral (effet pont ?)

Effets de contexte locaux - caractéristqiues intrinsèques / endogènes : normes locales partagées ex quartier pavillionnaire - caractéristiques exogènes : effet pont à proximité

Aller plus loin

GWR multiscalaire

Modèle mixte : GWR + modèle global dans le même modèle.

Nombre de variable déterminé pour chaque variable.

source("../R/gwr.multiscale_T.r") 

# On lance la MGWR
MGWR <- gwr.multiscale(formula = log_valeur_fonciere ~ log_dist_litt + 
                         log_surface_reelle_bati + log_surface_terrain,
                       data = dataSp, kernel = "exponential",
                       predictor.centered=rep(T, 3), # centrage des prédicteurs
                       adaptive = TRUE,
                       bws0 = rep(1,4)) # BW minimum pour l'optimisation
#> Error in is(data, "Spatial"): object 'dataSp' not found

mgwr.bw  <- round(MGWR[[2]]$bws,1) # Nombre de voisins pour chaque prédicteur
#> Error in eval(expr, envir, enclos): object 'MGWR' not found
mgwr.bw
#> Error in eval(expr, envir, enclos): object 'mgwr.bw' not found

# Exploration des résultats statistiques
print(MGWR)
#> Error in print(MGWR): object 'MGWR' not found

On constate que si l’effet de la surface bâtie sur les prix agit de façon assez globale à l’échelle de l’île, les deux autres prédicteurs agissent de manière beaucoup plus locales. Ainsi par exemple, l’effet de la distance au littoral relève d’un processus très localisé.

Cartographie des résultats

a <- mapGWR(MGWR$SDF, var = "log_dist_litt",var_TV = MGWR$SDF$log_dist_litt_TV,
       main.title = "Distance au littoral (bw = 77)")
#> Error in mapGWR(MGWR$SDF, var = "log_dist_litt", var_TV = MGWR$SDF$log_dist_litt_TV, : object 'MGWR' not found
b <- mapGWR(MGWR$SDF, var = "log_surface_reelle_bati",var_TV = MGWR$SDF$log_surface_reelle_bati_TV,
       main.title = "Surface bâtie (bw = 368)")
#> Error in mapGWR(MGWR$SDF, var = "log_surface_reelle_bati", var_TV = MGWR$SDF$log_surface_reelle_bati_TV, : object 'MGWR' not found
c <- mapGWR(MGWR$SDF, var = "log_surface_terrain",var_TV = MGWR$SDF$log_surface_terrain_TV,
       main.title = "Surface terrain (bw = 23)")
#> Error in mapGWR(MGWR$SDF, var = "log_surface_terrain", var_TV = MGWR$SDF$log_surface_terrain_TV, : object 'MGWR' not found

tmap_arrange(a,b,c)
#> Error in tmap_arrange(a, b, c): object 'a' not found

Régionalisation des sous-marchés immobiliers

On peut utiliser des coefficients locaux pour identifier des zones homogènes où on va essayer que les prédicteurs soient homogènes.

Sous-régions où les effets sont à peu près homogènes: minimisation de la variance.

Principe général: regrouper les observations qui se ressemblent (classification) mais aussi il faut qu’elles soient proches.

Plusieurs méthodes existent, nous allons utiliser SKATER qui utilise des graphes.

Méthode SKATER implémentée dans {spdep}.

Méthode en 4 étapes: - graphe de voisinage - pondération des arrètes (fonction de dissimilarité : GWR) - élagage de l’arbre - maximiser la variance inter-classe des graphes.

Première étape : préparation de la table des coefficients GWR

# Centrage-réduction pour rendre les coefficients comparables
gwrB.scaled <- GWR.exp$SDF %>% 
  as.data.frame() %>% 
  select(1:4) %>% 
  mutate_all(~scale(.)) %>% 
  rename_with(~paste(.x, "b", sep = "_"))
#> Error in rename_with(., ~paste(.x, "b", sep = "_")): could not find function "rename_with"

Deuxième étape : computation de l’algorithme SKATER

Définition du voisinage de chaque bien

knn <- knearneigh(GWR.exp$SDF, k=50) # 50 pour que le graphe soit fermé
#> Error in knearneigh(GWR.exp$SDF, k = 50): could not find function "knearneigh"
nb <- knn2nb(knn) # nombre de voisin
#> Error in knn2nb(knn): could not find function "knn2nb"
plot(nb, coords = coordinates(GWR.exp$SDF), col="blue")
#> Error in plot(nb, coords = coordinates(GWR.exp$SDF), col = "blue"): object 'nb' not found

Si graphe non fermé partout (connexité ?), il y a une erreur, augmentation du paramètre k.

Calibrage du coût des arêtes et de la pondération spatiale

costs <- nbcosts(nb, data = gwrB.scaled) # pondération des arêtes
#> Error in nbcosts(nb, data = gwrB.scaled): could not find function "nbcosts"
costsW <- nb2listw(nb, costs, style="B") # transformation en matrice de pondération, style : type de pondération, "B" binary
#> Error in nb2listw(nb, costs, style = "B"): could not find function "nb2listw"

Minimisation de l’arbre et classification

Suppression des arêtes qui ont des poids faibles (points qui se ressemblent pas). Au final, 1 point à une arête.

costsTree <- mstree(costsW)
#> Error in mstree(costsW): could not find function "mstree"
plot(costsTree, coords = coordinates(GWR.exp$SDF), col="blue", main = "Arbre portant minimal")
#> Error in plot(costsTree, coords = coordinates(GWR.exp$SDF), col = "blue", : object 'costsTree' not found

Découpe en 6 groupes

ncuts: nombre de groupe - 1

clus6 <- skater(edges = costsTree[,1:2], data = gwrB.scaled, ncuts = 5)
#> Error in skater(edges = costsTree[, 1:2], data = gwrB.scaled, ncuts = 5): could not find function "skater"

Troisième étape : analyse des résultats

Cartographie des clusters

dataClus <- dataSf %>% 
mutate(clus = as.factor(clus6$groups)) %>% 
bind_cols(gwrB.scaled)
#> Error in bind_cols(., gwrB.scaled): could not find function "bind_cols"

tmap_mode(mode = "view")
#> tmap mode set to interactive viewing
tm_shape(dataClus) + 
tm_symbols(col="clus", size=.8, palette = "Set1") +
tm_layout(title = "Classification en 6 groupes") 
#> Error in as.list.environment(environment()): object 'dataClus' not found

Caractérisation des clusters

nomVar <- c("log_dist_litt_b","log_surface_reelle_bati_b","log_surface_terrain_b","Intercept_b")

clusProfile <- dataClus[, c(nomVar, "clus")] %>% 
  group_by(clus) %>% 
  summarise_each(funs(mean)) %>% # moyenne des coeff GWR pour chaque cluster
  st_drop_geometry()
#> Error in summarise_each(., funs(mean)): could not find function "summarise_each"

clusLong <- reshape2::melt(clusProfile, id.vars = "clus")
#> Error in loadNamespace(x): there is no package called 'reshape2'

profilePlot <- ggplot(clusLong) +
  geom_bar(aes(x = variable, y = value), 
           fill = "grey25", 
           position = "identity", 
           stat = "identity") +
  scale_x_discrete("Effet") + 
  scale_y_continuous("Valeur moyenne par classe") +
  facet_wrap(~ clus) + 
  coord_flip() + 
  theme(strip.background = element_rect(fill="grey25"),
        strip.text = element_text(colour = "grey85", face = "bold"))
#> Error in ggplot(clusLong): could not find function "ggplot"

profilePlot
#> Error in eval(expr, envir, enclos): object 'profilePlot' not found

Atelier d’analyse spatiale

Thierry Feuillet