Skip to content

Commit

Permalink
Multiplecomparison temporary commit I should Have branched before.
Browse files Browse the repository at this point in the history
  • Loading branch information
MGousseff committed Oct 9, 2024
1 parent f15831e commit 86b8bd5
Show file tree
Hide file tree
Showing 6 changed files with 150 additions and 103 deletions.
16 changes: 16 additions & 0 deletions R/FactorialPlanWFs.R
Original file line number Diff line number Diff line change
@@ -0,0 +1,16 @@
library(FactoMineR)
library(Factoshiny)
# factorialPlanWfs<-function(intersec_sf, sfWfs = NULL){
# columnNames<-names(intersec_sf)
# LCZwfsNames<-grep( pattern = "LCZ*", x = names(intersec_sf), value = TRUE)


# }

test_sf_for_MCA<-st_drop_geometry(multicompare_test[
, grep( pattern = "LCZ*", x = names(multicompare_test), value = ) ])


test.MCA<-FactoMineR::MCA(X = test_sf_for_MCA, ncp = 4913)

Factoshiny::MCAshiny(test.MCA)
22 changes: 16 additions & 6 deletions R/LCZmodeByAgreementLevel.R
Original file line number Diff line number Diff line change
Expand Up @@ -11,16 +11,19 @@ LCZmodeByAgreementLevel <- function(intersec_sf, sfWfs = NULL){
LCZwfsNames<-grep( pattern = "LCZ*", x = names(intersec_sf), value = TRUE)
intersec_sf$LCZmode<-apply(intersec_sf[,LCZwfsNames], 1, Mode)

modeLCZSurfbyAgreement <- intersec_sf %>% group_by(maxAgree, LCZmode) %>% summarize(modeLCZsurf = sum(area)) %>% mutate(modeLCZSurfPerc = modeLCZsurf/sum(modeLCZsurf)*100)
modeLCZSurfbyAgreement <- intersec_sf %>% group_by(maxAgree, LCZmode) %>% summarize(modeLCZsurf = sum(drop_units(area))) %>% mutate(modeLCZSurfPerc = modeLCZsurf/sum(modeLCZsurf)*100)

generalProp<-intersec_sf %>%select(area, LCZmode) %>% mutate(totalArea=sum(area)) %>%
areaByAgreement <- intersec_sf %>% group_by(maxAgree) %>% summarize ( areaByAgree = sum(drop_units(area))) %>% mutate ( agreementSurfPerc = round(areaByAgree / sum(areaByAgree) * 100, digits = 3))

generalProp<-intersec_sf %>%select(area, LCZmode) %>% mutate(totalArea=sum(drop_units(area))) %>%
group_by(LCZmode) %>%
summarize(modeLCZGenSurfPerc = sum(area), totalArea = mean(totalArea)) %>%
mutate(modeLCZGenSurfPerc = modeLCZGenSurfPerc / totalArea *100 ) %>%
summarize(modeLCZGenSurfPerc = sum(drop_units(area)), totalArea = mean(totalArea)) %>%
mutate(modeLCZGenSurfPerc = round(modeLCZGenSurfPerc / totalArea *100 , digits = 3)) %>%
select(LCZmode, modeLCZGenSurfPerc)

modeLCZSurfbyAgreement<-left_join(modeLCZSurfbyAgreement, generalProp, by = "LCZmode") %>%
arrange(desc(maxAgree),desc(modeLCZSurfPerc))
modeLCZSurfbyAgreement<-left_join(modeLCZSurfbyAgreement, areaByAgreement, by = "maxAgree")

# if (!is.null(sfWfs)) {
# lengthSfWfs<-length(sfWfs)
Expand All @@ -36,7 +39,9 @@ LCZmodeByAgreementLevel <- function(intersec_sf, sfWfs = NULL){
# row.names(agreement_by_pair) <- compNames

# }
return(modeLCZSurfbyAgreement )


return(modeLCZSurfbyAgreement[,c("maxAgree", "agreementSurfPerc", "LCZmode", "modeLCZSurfPerc", "modeLCZGenSurfPerc", "modeLCZsurf", "areaByAgree")] )
}

Mode <- function(x) {
Expand All @@ -46,7 +51,12 @@ Mode <- function(x) {


LCZmodeTest <- LCZmodeByAgreementLevel(multicompare_test$intersec_sf)
LCZmodeTest[601:610,c(LCZwfsNamesTest,"LCZmode")]





LCZmodeTest[601:610,c(LCZwfsNamesTest,"LCZmode")]Je n'ai certes pas lu d'articles scientifiques à ce sujet, car je n'en ai pas le temps, ni sans doute la capacité
test1<-multicompare_test$intersec_sf
LCZwfsNamesTest<-grep( pattern = "LCZ*", x = names(test1), value = TRUE)
Expand Down
48 changes: 48 additions & 0 deletions R/agreementByModeLCZ.R
Original file line number Diff line number Diff line change
@@ -0,0 +1,48 @@
agreementByLCZlevel <- function(intersec_sf, sfWfs = NULL){
if ( !is.null(intersec_sf$geometry)) {
intersec_sf<-st_drop_geometry(intersec_sf)
}
columnNames<-names(intersec_sf)

LCZwfsNames<-grep( pattern = "LCZ*", x = names(intersec_sf), value = TRUE)
intersec_sf$LCZmode<-apply(intersec_sf[,LCZwfsNames], 1, Mode)

agreementByModeLCZ<- intersec_sf %>% group_by(LCZmode, maxAgree) %>% summarize(agreeByLCZsurf = sum(drop_units(area))) %>% mutate(LCZmode = as.integer(LCZmode), agreeByLCZPerc = round( agreeByLCZsurf/sum(agreeByLCZsurf)*100, digits = 2))

# areaByAgreement <- intersec_sf %>% group_by(maxAgree) %>% summarize ( areaByAgree = sum(drop_units(area))) %>% mutate ( agreementSurfPerc = round(areaByAgree / sum(areaByAgree) * 100, digits = 3))

generalProp<-intersec_sf %>% select(area, LCZmode) %>% mutate(totalArea=sum(drop_units(area))) %>%
group_by(LCZmode) %>%
summarize(modeLCZGenSurfPerc = sum(drop_units(area)), totalArea = mean(totalArea)) %>%
mutate(modeLCZGenSurfPerc = round(modeLCZGenSurfPerc / totalArea *100 , digits = 3)) %>%
select(LCZmode, modeLCZGenSurfPerc) %>% mutate(LCZmode = as.integer(LCZmode))

# modeLCZSurfbyAgreement<-left_join(modeLCZSurfbyAgreement, generalProp, by = "LCZmode") %>%
# arrange(desc(maxAgree),desc(modeLCZSurfPerc))
# modeLCZSurfbyAgreement<-left_join(modeLCZSurfbyAgreement, areaByAgreement, by = "maxAgree")

agreementByModeLCZ <- agreementByModeLCZ %>%
dplyr::arrange(LCZmode,desc(maxAgree)) %>% select(LCZmode, maxAgree, agreeByLCZPerc)

agreementByModeLCZ <- left_join(agreementByModeLCZ, generalProp, by = "LCZmode")

p<-ggplot(agreementByModeLCZ, aes(x = maxAgree, y = agreeByLCZPerc)) +
geom_bar(stat = "identity") +
facet_wrap(~ LCZmode)
print(p)
return(agreementByModeLCZ)
}

Mode <- function(x) {
ux <- unique(x)
unlist(ux[which.max(tabulate(match(x, ux)))])
}

library(units)
library(dplyr)
library(sf)
library(ggplot2)
library(lczexplore)
agreeByLCZmodeTest <- agreementByLCZlevel(multicompare_test$intersec_sf)

head(agreeByLCZmodeTest)
97 changes: 21 additions & 76 deletions R/compareMultipleLCZ.R
Original file line number Diff line number Diff line change
Expand Up @@ -20,99 +20,44 @@ compareMultipleLCZ<-function(sfList, LCZcolumns, refCrs=NULL, sfWf=NULL, trimPer
intersec_sf$area<-st_area(intersec_sf)
intersec_sf <- intersec_sf %>% subset(area>quantile(intersec_sf$area, probs=trimPerc) & !is.na(area))
print(nrow(intersec_sf))
intersec_sfnogeom<-st_drop_geometry(intersec_sf)
for (i in 1:(length(sfList) - 1)) {
for(j in (i+1):length(sfList)){
compName<-paste0(i,"_",j)
print(compName)
intersec_sfnogeom[,compName]<-intersec_sfnogeom[,i] == intersec_sfnogeom[,j]
}
}
rangeCol<-(length(sfList)+3):ncol(intersec_sfnogeom)
print(rangeCol)
# print(names(intersec_sfnogeom[,rangeCol]))
intersec_sfnogeom$nbAgree<-apply(intersec_sfnogeom[,rangeCol],MARGIN=1,sum)
intersec_sfnogeom$maxAgree<-apply(
X = intersec_sfnogeom[,1:length(sfList)], MARGIN = 1, function(x) max(table(x) ))
intersec_sf<-cbind(intersec_sfnogeom,intersec_sf$geometry) %>% st_as_sf()
intersec_sf
intersec_sfLong<-pivot_longer(intersec_sfnogeom,cols=rangeCol, names_to = "whichWfs", values_to = "agree")
intersec_sfLong$LCZref<-substr(intersec_sfLong$whichWfs,start = 1, stop=1 )
print(head(intersec_sfLong[,c(1,2,9:10)]))
whichLCZagree <- names(intersec_sfLong)[as.numeric(intersec_sfLong$LCZref)]
indRow<- seq_len(nrow(intersec_sfLong))
z<-data.frame(indRow, whichLCZagree)
intersec_sfLong$LCZvalue<-apply(z, 1, function(x) unlist(st_drop_geometry(intersec_sfLong)[x[1], x[2]]))
print(head(intersec_sfLong[,c(1,2,9:11)]))
nbWfs<-length(sfList)

intersec_sf<-computeAgreements(intersec_sf = intersec_sf, nbWfs = nbWfs)

return(intersec_sf)

output<-list(intersec_sf=intersec_sf, intersec_sfLong=intersec_sfLong)
}


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_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")
sf_IAU_auffargis <- importLCZvect("/home/gousseff/Documents/3_data/data_article_LCZ_diff_algos/IAU", file = "IAU_Auffargis.fgb", column = "lcz_primary")

sfList<-list(BDT11 = sfBDT_11_78030, BDT22 = sfBDT_22_78030, OSM11= sf_OSM_11_Auffargis, OSM22 = sf_OSM_22_Auffargis,
WUDAPT = sf_WUDAPT_78030)
sfList<-list(BDT11 = sfBDT_11_78030, BDT22 = sfBDT_22_78030, # OSM11= sf_OSM_11_Auffargis,
OSM22 = sf_OSM_22_Auffargis,
WUDAPT = sf_WUDAPT_78030, IAU = sf_IAU_auffargis)


multicompare_test<-compareMultipleLCZ(sfList = sfList, LCZcolumns = c(rep("LCZ_PRIMARY",4),"lcz_primary"),
sfWf = c("BDT11","BDT22","OSM11","OSM22","WUDAPT"),trimPerc = 0.5)
multicompare_test
multicompare_test<-compareMultipleLCZ(sfList = sfList, LCZcolumns = c(rep("LCZ_PRIMARY",3), rep("lcz_primary", 2)),
sfWf = c("BDT11","BDT22",
# "OSM11",
"OSM22","WUDAPT", "IAU"),trimPerc = 0.25)
multicompare_test %>% summary()


ggplot(data=multicompare_test$intersec_sf) +
require(ggplot2)
ggplot(data=multicompare_test) +
geom_sf(aes(fill=maxAgree, color=after_scale(fill)))+
scale_fill_gradient(low = "red" , high = "green", na.value = NA)

# hist(st_area(multicompare_test$intersec_sf$geometry))

# allSameList<-list(OSM11= sf_OSM_11_Auffargis, OSM11.2 = sf_OSM_11_Auffargis,
# OSM11.3 = sf_OSM_11_Auffargis, OSM11.4 = sf_OSM_11_Auffargis, OSM_11.5 = sf_OSM_11_Auffargis)
# showLCZ(sfList[[1]])

# sf_OSM_11_Auffargis[which.max(st_area(sf_OSM_11_Auffargis)),]

# max(st_area(sf_OSM_11_Auffargis))

# multicompare_test_all_same<-compareMultipleLCZ(sfList = allSameList,
# LCZcolumns = rep("LCZ_PRIMARY",5),
# sfWf = c("OSM1","OSM2", "OSM3", "OSM4", "OSM5"),
# trimPerc = 0.5)



# areas_test<-st_area(multicompare_test_all_same$intersec_sf)
# hist(areas_test)
# hist(st_area(sf_OSM_11_Auffargis$geometry))


# quantile(areas_test, prob = 0.5)
# test2<-multicompare_test_all_same$intersec_sf[
# st_area(multicompare_test_all_same$intersec_sf) ==
# max(st_area(multicompare_test_all_same$intersec_sf)),
# ]


# ggplot() +
# geom_sf(data=multicompare_test_all_same$intersec_sf, aes(fill=maxAgree))+
# scale_fill_gradient(low = "red" , high = "green", na.value = NA)

# ggplot() +
# geom_sf(data=test2, aes(color = "gray", fill=maxAgree)) +
# scale_fill_gradient(low = "red" , high = "green", na.value = NA) +
# scale_linewidth(range=c(8))
ggplot(data=multicompare_test) +
geom_sf(aes(fill=nbAgree, color=after_scale(fill)))+
scale_fill_gradient(low = "red" , high = "green", na.value = NA)

# ggplot() +
# geom_sf(data = sf_OSM_11_Auffargis[which.max(st_area(sf_OSM_11_Auffargis)),],
# aes(color = LCZ_PRIMARY)
# )
18 changes: 18 additions & 0 deletions R/computeAgreements.R
Original file line number Diff line number Diff line change
@@ -0,0 +1,18 @@
computeAgreements<-function(intersec_sf, nbWfs) {
intersec_sfnogeom<-st_drop_geometry(intersec_sf)
for (i in 1:(nbWfs - 1)) {
for(j in (i+1):nbWfs){
compName<-paste0(i,"_",j)
print(compName)
intersec_sfnogeom[,compName]<-intersec_sfnogeom[,i] == intersec_sfnogeom[,j]
}
} # TODO : refactor with a matrix reduction operation

rangeCol<-( (nbWfs+2):ncol(intersec_sfnogeom) ) # nbWfs+2 because of the area column
print(names(intersec_sfnogeom)[rangeCol])
# print(names(intersec_sfnogeom[,rangeCol]))
intersec_sfnogeom$nbAgree<-apply(intersec_sfnogeom[,rangeCol],MARGIN=1,sum)
intersec_sfnogeom$maxAgree<-apply(
X = intersec_sfnogeom[,1:nbWfs], MARGIN = 1, function(x) max(table(x) ))
intersec_sf<-cbind(intersec_sfnogeom,intersec_sf$geometry) %>% st_as_sf()
}
52 changes: 31 additions & 21 deletions R/createIntersect.R
Original file line number Diff line number Diff line change
Expand Up @@ -9,29 +9,39 @@ createIntersec<-function(sfList, LCZcolumns, refCrs=NULL, sfWf=NULL){
if (!is.null(sfWf) & length(sfWf) == length(sfList)){
names(intersec_sf)[1:(ncol(intersec_sf)-1)]<-paste0("LCZ",sfWf)
} else { names(intersec_sf)[1:(ncol(intersec_sf)-1)]<-paste0("LCZ",1:length(sfList)) }
intersec_sf

return(intersec_sf)

# intersec_sf[,1:(ncol(intersec_sf)-1)]
# test<- apply(
# X = st_drop_geometry(intersec_sf)[,1:(ncol(intersec_sf)-1)],
# MARGIN = 2, FUN = as.factor)
# return(test)

}

# sfBDT_11_78030<-importLCZvect(dirPath="/home/gousseff/Documents/0_DocBiblioTutosPublis/0_ArticlesScientEtThèses/ArticleComparaisonLCZGCWUDAPTEXPERTS/BDT/2011/bdtopo_2_78030",
# file="rsu_lcz.fgb", column="LCZ_PRIMARY")
# class(sfBDT_11_78030)
# sfBDT_22_78030<-importLCZvect(dirPath="/home/gousseff/Documents/0_DocBiblioTutosPublis/0_ArticlesScientEtThèses/ArticleComparaisonLCZGCWUDAPTEXPERTS/BDT/2022/bdtopo_3_78030",
# file="rsu_lcz.fgb", column="LCZ_PRIMARY")
# sf_OSM_11_Auffargis<-importLCZvect(dirPath="/home/gousseff/Documents/0_DocBiblioTutosPublis/0_ArticlesScientEtThèses/ArticleComparaisonLCZGCWUDAPTEXPERTS/OSM/2011/osm_Auffargis/",
# file="rsu_lcz.fgb", column="LCZ_PRIMARY")
# sf_OSM_22_Auffargis<-importLCZvect(dirPath="/home/gousseff/Documents/0_DocBiblioTutosPublis/0_ArticlesScientEtThèses/ArticleComparaisonLCZGCWUDAPTEXPERTS/OSM/2022/osm_Auffargis/",
# file="rsu_lcz.fgb", column="LCZ_PRIMARY")
# sf_WUDAPT_78030<-importLCZvect("/home/gousseff/Documents/0_DocBiblioTutosPublis/0_ArticlesScientEtThèses/ArticleComparaisonLCZGCWUDAPTEXPERTS/WUDAPT",
# file ="wudapt_78030.geojson", column="lcz_primary")
#
# 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<-createIntersec(sfList = sfList, LCZcolumns = c(rep("LCZ_PRIMARY",4),"lcz_primary"),
# sfWf = c("BDT11","BDT22","OSM11","OSM22","WUDAPT"))
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")
sf_IAU_auffargis <- importLCZvect("/home/gousseff/Documents/3_data/data_article_LCZ_diff_algos/IAU", file = "IAU_Auffargis.fgb", column = "lcz_primary")

sfList<-list(BDT11 = sfBDT_11_78030, BDT22 = sfBDT_22_78030, OSM11= sf_OSM_11_Auffargis, OSM22 = sf_OSM_22_Auffargis,
WUDAPT = sf_WUDAPT_78030, IAU = sf_IAU_auffargis)



intersected<-createIntersec(sfList = sfList, LCZcolumns = c(rep("LCZ_PRIMARY",4),rep("lcz_primary",2)),
sfWf = c("BDT11","BDT22","OSM11","OSM22","WUDAPT", "IAU"))
is.factor(sf_OSM_22_Auffargis$LCZ_PRIMARY)

#
#
# test_list<-list(a=c(1,2),b="top",c=TRUE)
Expand Down

0 comments on commit 86b8bd5

Please sign in to comment.