diff --git a/R/compareLCZ.R b/R/compareLCZ.R index 8d949f6..7b0a4f8 100644 --- a/R/compareLCZ.R +++ b/R/compareLCZ.R @@ -365,7 +365,7 @@ compareLCZ<-function(sf1,geomID1="",column1,confid1="",wf1="bdtopo_2_2", write.table(x=echIntExpo, file =nom, append = TRUE, quote = TRUE, sep = ";", eol = "\n", na = "NA", dec = ".", qmethod = c("escape", "double"), - fileEncoding = "", col.names=FALSE,row.names=F) + fileEncoding = "", col.names=FALSE, row.names=F) } } ################################################### @@ -393,12 +393,6 @@ areas<-matConfOut$areas percAgg<-matConfOut$percAgg - - - - - - ################################################ # GRAPHICS ################################################ diff --git a/R/compareMultipleLCZ.R b/R/compareMultipleLCZ.R index 7225e26..f796793 100644 --- a/R/compareMultipleLCZ.R +++ b/R/compareMultipleLCZ.R @@ -21,9 +21,7 @@ compareMultipleLCZ<-function(sfList, LCZcolumns, refCrs=NULL, sfWf=NULL, trimPerc=0.05){ echInt<-createIntersect(sfList = sfList, columns = LCZcolumns , refCrs= refCrs, sfWf = sfWf) print(nrow(echInt)) - echInt$area<-st_area(echInt) echInt <- echInt %>% subset(area>quantile(echInt$area, probs=trimPerc) & !is.na(area)) - print(nrow(echInt)) echIntnogeom<-st_drop_geometry(echInt) for (i in 1:(length(sfList) - 1)) { for(j in (i+1):length(sfList)){ diff --git a/R/createIntersect.R b/R/createIntersect.R index 5e64dbf..a4e00a0 100644 --- a/R/createIntersect.R +++ b/R/createIntersect.R @@ -10,7 +10,7 @@ #' are assigned to geometries resulting from intersection of all input geometries #' @export #' @examples -createIntersect<-function(sfList, columns, refCrs=NULL, sfWf=NULL){ +createIntersect<-function(sfList, columns, refCrs=NULL, sfWf=NULL, minZeroArea=0.0001){ sfInt<-sfList[[1]] %>% select(columns[1]) if (is.null(refCrs)){refCrs<-st_crs(sfInt)} for (i in 2:length(sfList)){ @@ -19,8 +19,10 @@ createIntersect<-function(sfList, columns, refCrs=NULL, sfWf=NULL){ sfInt<-st_intersection(sfInt,sfProv) } if (!is.null(sfWf) & length(sfWf) == length(sfList)){ - names(sfInt)[1:(ncol(sfInt)-1)]<-paste0("LCZ",sfWf) + names(sfInt)[1:(ncol(sfInt)-1)]<-sfWf } else { names(sfInt)[1:(ncol(sfInt)-1)]<-paste0("LCZ",1:length(sfList)) } + sfInt$area<-units::drop_units(st_area(sfInt$geometry)) + sfInt<-sfInt[sfInt$area>minZeroArea,] return(sfInt) } diff --git a/R/multipleCramer.R b/R/multipleCramer.R index fd3a14a..fdc4d5c 100644 --- a/R/multipleCramer.R +++ b/R/multipleCramer.R @@ -1,12 +1,15 @@ -#' Computes Cramer's V for all pairs of one hot encoded levels of input LCZ classifications -#' @param sfInt an sf object with several LCZ classifications to compare on the same intersected geometries, -#' typically an output of the createIntersec function +#' Computes Cramer's V for all pairs of onehot encoded levels of input LCZ classifications +#' @param sfInt an sf object with several LCZ classifications on the same (intersected) geometries, +#' typically an output of the createIntersect function #' @param columns a vector which contains names of LCZ classification columns +#' @param nbOutAssociations the number of significant associations we want to extract, ie pair of levels +#' from two different workflows whose Cramer's V is high (association between LCZ 101 from workflow 1 +#' and LCZ 101 from workflow 2 are exclude, for instance, but can be read from the cramerLong output) #' @importFrom ggplot2 geom_sf guides ggtitle aes #' @importFrom caret dummyVars #' @import sf dplyr cowplot forcats units tidyr RColorBrewer utils grDevices rlang -#' @return an sf file with values of LCZ from all the input -#' are assigned to geometries resulting from intersection of all input geometries +#' @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 multipleCramer<-function(sfInt, columns, nbOutAssociations ){ @@ -27,29 +30,37 @@ multipleCramer<-function(sfInt, columns, nbOutAssociations ){ # print(table(sfDummy[,c(i,j)])) } } + rownames(Vs)<-names(sfDummy) colnames(Vs)<-names(sfDummy) - threshold <- Vs %>% as.vector %>% unique %>% sort(decreasing=TRUE) %>% head(n = (nbOutAssociations+1)) %>% min - print(paste0("threshold = ", threshold)) + VsLong<-data.frame( + LCZ_type_1 = rownames(Vs), + LCZ_type_2 = rep(colnames(Vs), each = nrow(Vs)), + cramerVs = as.vector(Vs)) %>% arrange(desc(cramerVs))%>% na.omit() - Mask<-Vs - Mask<-(Mask>=threshold) - Mask[Vs0.004 & testMask<1) -testMask[!testMask]<-NA - -keptRows<-apply( - X = apply(X = testMask, MARGIN = 1, is.na), - MARGIN = 1, sum)=1 | signifAssoc<0.004)]<-NA +testVs$cramerLong %>% head(10) + +intersectedDf<-st_drop_geometry(intersected) +str(intersectedDf) +summary(intersectedDf$area) +min(intersectedDf$area) +dataTest<-intersectedDf[ ,-6] +dataTest <- as.data.frame(lapply(X = dataTest, factor)) +summary(dataTest) +str(dataTest) +weights<-(intersectedDf$area/sum(intersectedDf$area)) +length(weights) +auffargisMCA<-MCA(X = dataTest[, names(dataTest)!="area"], ncp = 5, row.w = weights) +plot.MCA(auffargisMCA, invisible = c("ind")) + +dataTestNo107<-dataTest[apply(dataTest, 1, function(x) all(x!="107")),] +nrow(dataTestNo107) +weightsNo107<-(intersectedDf$area/sum(intersectedDf$area))[ + apply(dataTest, 1, function(x) all(x!="107"))] +length(weightsNo107) + +auffargisMCANo107<-MCA(X = dataTestNo107[,names(dataTest)!="area"], ncp = 5) #, row.w = weightsNo107) +plot.MCA(auffargisMCANo107, invisible= c("ind")) +auffargisMCANo107Weights<-MCA(X = dataTestNo107[,names(dataTest)!="area"], ncp = 5, row.w = weightsNo107) +plot.MCA(auffargisMCANo107Weights, invisible= c("ind")) + + +str(auffargisMCANo107) +