Skip to content

Commit

Permalink
new functions for LCZ compare article
Browse files Browse the repository at this point in the history
  • Loading branch information
MGousseff committed Dec 30, 2024
1 parent 92e1552 commit 36557f2
Show file tree
Hide file tree
Showing 6 changed files with 141 additions and 52 deletions.
44 changes: 44 additions & 0 deletions R/ConcatAllLocationsAllWf.R
Original file line number Diff line number Diff line change
@@ -0,0 +1,44 @@
#' In a given list of directories the function looks for LCZ datafiles, return a datasets LCZ values for each geom and workflow
#' @param dirList the list of directories for which the different LCZ files will be intersected
#' @param workflowNames sets the names of workflows and define the name of the files which will be loaded and intersected
#' @param for each diretory from dirList, a location name must be fed to the function
#' @importFrom ggplot2 geom_sf guides ggtitle aes
#' @import sf dplyr cowplot forcats units tidyr RColorBrewer utils grDevices rlang
#' @return returns graphics of comparison and an object called matConfOut which contains :
#' matConfLong, a confusion matrix in a longer form,
#' matConfPlot is a ggplot2 object showing the confusion matrix.
#' percAgg is the general agreement between the two sets of LCZ, expressed as a percentage of the total area of the study zone
#' pseudoK is a heuristic estimate of a Cohen's kappa coefficient of agreement between classifications
#' If saveG is not an empty string, graphics are saved under "saveG.png"
#' @export
#' @examples
#'
concatIntersectedLocations<-function(dirList, locations, workflowNames = c("osm","bdt","iau","wudapt")){
allLocAllWfSf<-matrix(ncol = 5, nrow = 0)
allLocAllWfSf<-as.data.frame(allLocAllWfSf)
names(allLocAllWfSf)<- c("lcz_primary", "location", "wf", "area", "geometry")
for( i in 1:length(dirList)){
dirPath<-dirList[1]
aLocation<-locations[i]
sfList<-loadMultipleSf(dirPath = dirPath,
workflowNames = workflowNames , location = aLocation )
if(substr(dirPath, nchar(dirPath), nchar(dirPath))!="/"){dirPath<-paste0(dirPath, "/")}
zoneSfPath<-paste0(dirPath,"zone.fgb")
zoneSf<-read_sf(zoneSfPath)
sfList<-repairRoadsIAU(sfList = sfList, zoneSf = zoneSf, location = location)
concatSf<-concatAlocationWorkflows(sfList = sfList,
location = location, refCrs = 1)
allLocAllWfSf<-rbind(allLocAllWfSf, concatSf)
}

allLocAllWfSf<- st_as_sf(allLocAllWfSf)
}

rootDir<-"/home/gousseff/Documents/3_data/data_article_LCZ_diff_algos/newDataTree"
setwd(rootDir)
allLCZDirNames <- list.dirs()[-1]
allLCZDirNames <- substr(allLCZDirNames, start = 2, stop = 1000)
allLocationsNames<-substr(allLCZDirNames, start = 2, stop = 1000)
allLCZDirNames<-paste0(rootDir, allLCZDirNames, "/")
allLocallWfs<-concatIntersectedLocations(
dirList = allLCZDirNames, locations = allLocationsNames , workflowNames = c("osm","bdt","iau","wudapt"))
4 changes: 2 additions & 2 deletions R/compareMultipleLCZ.R
Original file line number Diff line number Diff line change
Expand Up @@ -49,9 +49,9 @@ compareMultipleLCZ<-function(sfInt, LCZcolumns, sfWf=NULL, trimPerc=0.05){
whichLCZagree <- gsub( x = sfIntLong$whichWfs, pattern = "(.*)(_)(.*)", replacement = "\\1")
indRow<- seq_len(nrow(sfIntLong))
z<-data.frame(indRow, whichLCZagree)
print(head(z))

sfIntLong$LCZvalue<-apply(z, 1, function(x) unlist(st_drop_geometry(sfIntLong)[x[1], x[2]]))
print(head(sfIntLong[,c(1,2,9:11)]))

sfInt<-cbind(sfIntnogeom,sfInt$geometry) %>% st_as_sf()


Expand Down
5 changes: 4 additions & 1 deletion R/importLCZvect.R
Original file line number Diff line number Diff line change
Expand Up @@ -48,7 +48,10 @@ importLCZvectFromFile<-function(dirPath, file="rsu_lcz.geojson", column, geomID=
colErr<-c("It seems that some of the columns you try to import do not exist in the source file,
are you sure you meant ",
paste(badCol),"?")
if (prod(inCol)==0){ stop(colErr) } else {
if (prod(inCol)==0){


stop(colErr) } else {
if (drop== TRUE) {sfFile<-sf::st_read(dsn=fileName,quiet=!verbose)[,colonnes] } else {
sfFile<-sf::st_read(dsn=fileName,quiet=!verbose)[,]}
}
Expand Down
1 change: 1 addition & 0 deletions R/loadMultipleSf.R
Original file line number Diff line number Diff line change
Expand Up @@ -19,6 +19,7 @@ loadMultipleSf<-function(dirPath, workflowNames = c("osm","bdt","iau","wudapt"),
"101"="101","102"="102","103"="103","104"="104", "105"="105","106"="106","107"="107",
"101"="11","102"="12","103"="13","104"="14", "105"="15", "106"="16","107"="17",
"101"="A","102"="B","103"="C","104"="D","105"="E","106"="F","107"="G")

sfList<-list()
for (i in workflowNames){
inName<-paste0(dirPath, i, "_lcz.fgb")
Expand Down
52 changes: 52 additions & 0 deletions R/plotLCZaLocation.R
Original file line number Diff line number Diff line change
@@ -0,0 +1,52 @@
#' in a given directory, if several LCZ files are present, it plots the repartition of LCZ regarding of their source (workflow)
#' NOTE: to represent the map of LCZ for a given fil, use `showLCZ` function instead
#' @param workflownames is a vector of prefixes. The LCZ files must benames workflow_rsu.fgb
#' @param plotNow. If TRUE, the boxplot of the repartition will be printed
#' @param plotSave. If TRUE, the plot will be saved in the directory pointed by dirPath
#' @importFrom ggplot2 geom_sf guides ggtitle aes
#' @importFrom caret dummyVars
#' @import sf dplyr cowplot forcats units tidyr RColorBrewer utils grDevices rlang
#' @return Cramer's V between pairs of levels, in a matrix (cramerMatrix) or long form (cramerLong),
#' and a dataframe with the nbOutAssociation most significant association
#' @export
#' @examples
plotLCZaLocation<-function(dirPath, location, workflowNames = c("osm","bdt","iau","wudapt"),
plotNow = FALSE, plotSave = TRUE){
colorMap<-c("#8b0101","#cc0200","#fc0001","#be4c03","#ff6602","#ff9856",
"#fbed08","#bcbcba","#ffcca7","#57555a","#006700","#05aa05",
"#648423","#bbdb7a","#010101","#fdf6ae","#6d67fd")
names(colorMap)<-c(1:10,101:107)
etiquettes<-c("LCZ 1: Compact high-rise","LCZ 2: Compact mid-rise","LCZ 3: Compact low-rise",
"LCZ 4: Open high-rise","LCZ 5: Open mid-rise","LCZ 6: Open low-rise",
"LCZ 7: Lightweight low-rise","LCZ 8: Large low-rise",
"LCZ 9: Sparsely built","LCZ 10: Heavy industry",
"LCZ A: Dense trees", "LCZ B: Scattered trees",
"LCZ C: Bush,scrub","LCZ D: Low plants",
"LCZ E: Bare rock or paved","LCZ F: Bare soil or sand","LCZ G: Water")

sfList<-loadMultipleSf(dirPath = dirPath,
workflowNames = workflowNames , location = location )
if(substr(dirPath, nchar(dirPath), nchar(dirPath))!="/"){dirPath<-paste0(dirPath, "/")}
zoneSfPath<-paste0(dirPath,"zone.fgb")
zoneSf<-read_sf(zoneSfPath)
sfList<-repairRoadsIAU(sfList = sfList, zoneSf = zoneSf, location = location)
concatSf<-concatAlocationWorkflows(sfList = sfList,
location = location, refCrs = 1)

surfaces<-concatSf %>%
mutate(wf = factor(wf, levels = c("bdt", "osm", "wudapt", "iau"))) %>%
group_by(wf, lcz_primary) %>% summarise(area=sum(area), location=unique(location))

location<-unique(surfaces$location)
outPlot<-ggplot(surfaces, aes(fill=lcz_primary, y=area, x=wf)) +
geom_col() +
# scale_fill_viridis(discrete = T) +
scale_fill_manual(values=colorMap, breaks = names(colorMap), labels = etiquettes) +
ggtitle(paste0("LCZ repartition by workflow for ", location))

plotName<-paste0(dirPath, "LCZbyWfBoxplot.png")
if(plotSave){ggsave(plotName, outPlot)}
if (plotNow){print(outPlot)}

return(outPlot)
}
87 changes: 38 additions & 49 deletions inst/tinytest/test_compareMultiple.R
Original file line number Diff line number Diff line change
Expand Up @@ -7,60 +7,49 @@ library(ggplot2)
library(cowplot)
library(forcats)

sfBDT_11_78030<-importLCZvect(dirPath="/home/gousseff/Documents/3_data/data_article_LCZ_diff_algos/GeoClimate/2011/bdtopo_2_78030/",
file="rsu_lcz.fgb", column="LCZ_PRIMARY")
class(sfBDT_11_78030)
sfBDT_22_78030<-importLCZvect(dirPath="/home/gousseff/Documents/3_data/data_article_LCZ_diff_algos/GeoClimate/2022/bdtopo_3_78030/",
file="rsu_lcz.fgb", column="LCZ_PRIMARY")
sf_OSM_11_Auffargis<-importLCZvect(dirPath="//home/gousseff/Documents/3_data/data_article_LCZ_diff_algos/GeoClimate/2011/osm_Auffargis/",
file="rsu_lcz.fgb", column="LCZ_PRIMARY")
sf_OSM_22_Auffargis<-importLCZvect(dirPath="/home/gousseff/Documents/3_data/data_article_LCZ_diff_algos/GeoClimate/2022/osm_Auffargis/",
file="rsu_lcz.fgb", column="LCZ_PRIMARY")
sf_WUDAPT_78030<-importLCZvect("/home/gousseff/Documents/3_data/data_article_LCZ_diff_algos/WUDAPT/",
file ="wudapt_Auffargis.fgb", column="lcz_primary")
sfList<-loadMultipleSf(dirPath = "/home/gousseff/Documents/3_data/data_article_LCZ_diff_algos/newDataTree/Drancy/",
workflowNames = c("osm","bdt","iau","wudapt"), location = "Drancy" )

sfList<-list(BDT11 = sfBDT_11_78030, BDT22 = sfBDT_22_78030, OSM11= sf_OSM_11_Auffargis, OSM22 = sf_OSM_22_Auffargis,
WUDAPT = sf_WUDAPT_78030)
# showLCZ(sfList[[1]])



intersected<-createIntersect(sfList = sfList, columns = c(rep("LCZ_PRIMARY",4),"lcz_primary"),
sfWf = c("BDT11","BDT22","OSM11","OSM22","WUDAPT"))
names(intersected) %>% reformulate()


# test_list<-list(a=c(1,2),b="top",c=TRUE)
# length(test_list)
# for (i in test_list[2:3]) print(str(i))
intersected<-createIntersect(sfList = sfList, columns = rep("lcz_primary", 4),
sfWf = c("osm","bdt","iau","wudapt"))

multicompare_test<-compareMultipleLCZ(intersected,
LCZcolumns = c("BDT11","BDT22","OSM11","OSM22","WUDAPT"),trimPerc = 0.5)
LCZcolumns = c("osm","bdt","iau","wudapt"),
trimPerc = 0.5)

testAreas<-workflowAgreeAreas(multicompare_test$sfIntLong)
testAreas$disagreeAreas
testAreas$agreeAreas

multicompare_test

test<-multicompare_test$sfIntLong
test2<-test %>% subset(agree==TRUE) %>% group_by(LCZvalue) %>% summarize(agreementArea=sum(area)) %>%
mutate(percAgreementArea=agreementArea/sum(agreementArea))

testWfAgree<-test %>% subset(agree==TRUE) %>% group_by(whichWfs) %>% summarize(agreementArea=sum(area))

test<-multicompare_test$sfInt[,1:5] %>% st_drop_geometry()
prov1<-apply(X = test, MARGIN = 1, table )
prov2<-apply(X = test, MARGIN = 1, function(x) max(table(x)) )

head(prov1)
head(prov2)

plot1<-showLCZ(sf = multicompare_test$sfInt, column="BDT22", wf="BDT22")
plot2<-showLCZ(sf = multicompare_test$sfInt, column="BDT11", wf="BDT1111")
plot3<-showLCZ(sf = multicompare_test$sfInt, column="OSM22", wf="OSM22")
plot4<-showLCZ(sf = multicompare_test$sfInt, column="WUDAPT", wf="WUDAPT")
plot5<-ggplot(data=multicompare_test$sfInt) +
geom_sf(aes(fill=nbAgree, color=after_scale(fill)))+
scale_fill_gradient(low = "red" , high = "green", na.value = NA)
plot_grid(plot1, plot2, plot3, plot4, plot5)
osm<-importLCZvect(dirPath = "/home/gousseff/Documents/3_data/data_article_LCZ_diff_algos/newDataTree/Dourdan/",
file = "osm_lcz.fgb")
bdt<-importLCZvect(dirPath = "/home/gousseff/Documents/3_data/data_article_LCZ_diff_algos/newDataTree/Dourdan/",
file = "bdt_lcz.fgb")
iau<-importLCZvect(dirPath = "/home/gousseff/Documents/3_data/data_article_LCZ_diff_algos/newDataTree/Dourdan/",
file = "iau_lcz.fgb")
wudapt<-importLCZvect(dirPath = "/home/gousseff/Documents/3_data/data_article_LCZ_diff_algos/newDataTree/Dourdan/",
file = "wudapt_lcz.fgb")

# multicompare_test
#
# test<-multicompare_test$sfIntLong
# test2<-test %>% subset(agree==TRUE) %>% group_by(LCZvalue) %>% summarize(agreementArea=sum(area)) %>%
# mutate(percAgreementArea=agreementArea/sum(agreementArea))
#
# testWfAgree<-test %>% subset(agree==TRUE) %>% group_by(whichWfs) %>% summarize(agreementArea=sum(area))
#
# test<-multicompare_test$sfInt[,1:5] %>% st_drop_geometry()
# prov1<-apply(X = test, MARGIN = 1, table )
# prov2<-apply(X = test, MARGIN = 1, function(x) max(table(x)) )
#
# head(prov1)
# head(prov2)
#
# plot1<-showLCZ(sf = multicompare_test$sfInt, column="BDT22", wf="BDT22")
# plot2<-showLCZ(sf = multicompare_test$sfInt, column="BDT11", wf="BDT1111")
# plot3<-showLCZ(sf = multicompare_test$sfInt, column="OSM22", wf="OSM22")
# plot4<-showLCZ(sf = multicompare_test$sfInt, column="WUDAPT", wf="WUDAPT")
# plot5<-ggplot(data=multicompare_test$sfInt) +
# geom_sf(aes(fill=nbAgree, color=after_scale(fill)))+
# scale_fill_gradient(low = "red" , high = "green", na.value = NA)
# plot_grid(plot1, plot2, plot3, plot4, plot5)

0 comments on commit 36557f2

Please # to comment.