diff --git a/R/FactorialPlanWFs.R b/R/FactorialPlanWFs.R new file mode 100644 index 0000000..ac14cba --- /dev/null +++ b/R/FactorialPlanWFs.R @@ -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) diff --git a/R/LCZmodeByAgreementLevel.R b/R/LCZmodeByAgreementLevel.R index c743ef4..c9fe64c 100644 --- a/R/LCZmodeByAgreementLevel.R +++ b/R/LCZmodeByAgreementLevel.R @@ -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) @@ -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) { @@ -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) diff --git a/R/agreementByModeLCZ.R b/R/agreementByModeLCZ.R new file mode 100644 index 0000000..4d1d03f --- /dev/null +++ b/R/agreementByModeLCZ.R @@ -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) diff --git a/R/compareMultipleLCZ.R b/R/compareMultipleLCZ.R index 3682538..f1f3c02 100644 --- a/R/compareMultipleLCZ.R +++ b/R/compareMultipleLCZ.R @@ -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) -# ) diff --git a/R/computeAgreements.R b/R/computeAgreements.R new file mode 100644 index 0000000..3442331 --- /dev/null +++ b/R/computeAgreements.R @@ -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() +} diff --git a/R/createIntersect.R b/R/createIntersect.R index 9f26180..73a3217 100644 --- a/R/createIntersect.R +++ b/R/createIntersect.R @@ -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)