############################################################
#Downloading appropriate libraries if not already downloaded
############################################################
list.of.packages <- c("landscapemetrics","raster","RCurl")
new.packages <- list.of.packages[!(list.of.packages %in% installed.packages()[,"Package"])]
if(length(new.packages)) install.packages(new.packages)
install.packages("landscapemetrics")
## package 'landscapemetrics' successfully unpacked and MD5 sums checked
##
## The downloaded binary packages are in
## C:\Users\homatt11\AppData\Local\Temp\RtmpGmeDP4\downloaded_packages
library(RCurl)
install.packages("shiny")
## package 'shiny' successfully unpacked and MD5 sums checked
##
## The downloaded binary packages are in
## C:\Users\homatt11\AppData\Local\Temp\RtmpGmeDP4\downloaded_packages
library(shiny)
############################################################
############################################################
require(landscapemetrics)
## Loading required package: landscapemetrics
require(raster)
## Loading required package: raster
## Loading required package: sp
require(RCurl)
############################################################
############################################################
#
# Read forest data
#
############################################################
forest=as.matrix(read.table("forest_1x1.txt"))
forest.r=raster(forest) #georeference the data
#Plotting forest
image(forest.r, main="Forest", col=grey(1:0))
#Calculating the window size (Size of graph)
forest.lengthOfArea = ncol(forest)
forest.widthOfArea = nrow(forest)
############################################################
# Forest output
############################################################
forest.total.area = lsm_c_ca(forest.r, directions = 8)$value[2]*forest.lengthOfArea*forest.widthOfArea*10000 # eight neightbour rule
forest.prop.area = lsm_c_ca(forest.r, directions = 8)$value[2]*10000
forest.number.patches = lsm_c_np(forest.r, directions = 8)$value[2]
forest.total.edge = lsm_l_te(forest.r, count_boundary = FALSE)$value/(10000)*forest.lengthOfArea*forest.widthOfArea*31.0461
forest.total.core.area = lsm_c_tca(forest.r, directions = 8, consider_boundary = FALSE, edge_depth = 1)$value[2]*forest.lengthOfArea*forest.widthOfArea*10000
forest.prop.core.area = lsm_c_tca(forest.r, directions = 8, consider_boundary = FALSE, edge_depth = 1)$value[2]*10000
forest.intermediate <- lsm_p_area(forest.r, directions = 8)
forest.max.patch.area <- max(forest.intermediate[which(forest.intermediate$class==1),]$value)*
forest.lengthOfArea*forest.widthOfArea*10000
forest.max.patch.index <- max(forest.intermediate[which(forest.intermediate$class==1),]$value)*10000
forest.mean.patch.area <- forest.total.area/forest.number.patches
forest.totalArea <- paste("Forest total area: ", forest.total.area)
forest.propTotalArea <- paste("Forest proportion total area: ", round(forest.prop.area,4))
forest.numberPatches <- paste("Forest number of patches: ", forest.number.patches)
forest.totalEdge <- paste("Forest total edge: ", round(forest.total.edge))
forest.totalCoreArea <- paste("Forest total core area: ", forest.total.core.area)
forest.propCoreArea <-paste("Forest proportion core area: ", round(forest.prop.core.area,4))
forest.maxPatchArea <- paste("Forest max patch area: ", forest.max.patch.area)
forest.maxPatchIndex <- paste("Forest max patch index: ", round(forest.max.patch.index,4))
forest.meanPatchArea <- paste("Forest mean patch area: ", round(forest.mean.patch.area,4))
forest.metrics=c(forest.totalArea,forest.propTotalArea,forest.numberPatches,
forest.totalEdge,forest.totalCoreArea,forest.propCoreArea,
forest.maxPatchArea,forest.maxPatchIndex,forest.meanPatchArea)
forest.metrics
## [1] "Forest total area: 21500"
## [2] "Forest proportion total area: 0.2123"
## [3] "Forest number of patches: 66"
## [4] "Forest total edge: 7720"
## [5] "Forest total core area: 14590"
## [6] "Forest proportion core area: 0.1441"
## [7] "Forest max patch area: 8225"
## [8] "Forest max patch index: 0.0812"
## [9] "Forest mean patch area: 325.7576"
############################################################
# For forest versus core area
############################################################
f.intermediate <- lsm_p_area(forest.r, directions = 8)
forest.patch.area <- f.intermediate[which(f.intermediate$class==1),]$value*forest.lengthOfArea*forest.widthOfArea*10000
c.intermediate <- lsm_p_core(forest.r, directions = 8, consider_boundary = FALSE, edge_depth = 1)
forest.path.core.area <- c.intermediate[which(c.intermediate$class==1),]$value*forest.lengthOfArea*forest.widthOfArea*10000
plot(forest.patch.area,forest.path.core.area,main="Forest: patch size area vs core area")
test<-cor.test(forest.patch.area,forest.path.core.area)
testMethod <- paste("Forest versus Core Area:", test$method)
testAlternative <- paste("Alternative hypothesis: ", test$alternative)
tDf <- paste("Degrees of freedom: ", test$parameter)
tValue <- paste("t-value: ", round(test$statistic,4))
tpValue <- paste("p-value: ", round(test$p.value,4))
correlationValue <- paste("Correlation: ", round(test$estimate,4))
forest.core.correlation=c(testMethod,testAlternative,tDf,
tValue,tpValue,correlationValue)
forest.core.correlation
## [1] "Forest versus Core Area: Pearson's product-moment correlation"
## [2] "Alternative hypothesis: two.sided"
## [3] "Degrees of freedom: 64"
## [4] "t-value: 444.2706"
## [5] "p-value: 0"
## [6] "Correlation: 0.9998"
############################################################
# For forest versus perimeter
############################################################
p.intermediate <- lsm_p_perim(forest.r, directions = 8)
forest.patch.perim <- p.intermediate[which(p.intermediate$class==1),]$value*312.5
plot(forest.patch.area,forest.patch.perim,main="Forest: patch size area vs perimeter")
test<-cor.test(forest.patch.area,forest.patch.perim)
testMethod <- paste("Forest versus Perimeter Correlation:", test$method)
testAlternative <- paste("Alternative hypothesis: ", test$alternative)
tDf <- paste("Degrees of freedom: ", test$parameter)
tValue <- paste("t-value: ", round(test$statistic,4))
tpValue <- paste("p-value: ", round(test$p.value,4))
correlationValue <- paste("Correlation: ", round(test$estimate,4))
forest.perimeter=c(testMethod,testAlternative,tDf,
tValue,tpValue,correlationValue)
forest.perimeter
## [1] "Forest versus Perimeter Correlation: Pearson's product-moment correlation"
## [2] "Alternative hypothesis: two.sided"
## [3] "Degrees of freedom: 64"
## [4] "t-value: 143.2683"
## [5] "p-value: 0"
## [6] "Correlation: 0.9984"
############################################################
# For forest versus perimeter area ratio
############################################################
intermediate <- lsm_p_para(forest.r, directions = 8)
forest.patch.perim.area.ratio <- intermediate[which(intermediate$class==1),]$value
plot(forest.patch.area,forest.patch.perim.area.ratio,main="Forest: patch size area vs perimeter area ratio")
test<-cor.test(forest.patch.area,forest.patch.perim.area.ratio,method="spearman")
## Warning in cor.test.default(forest.patch.area, forest.patch.perim.area.ratio, :
## Cannot compute exact p-value with ties
testMethod <- paste("Forest versus Perimeter/Area: ", test$method)
testAlternative <- paste("Alternative hypothesis: ", test$alternative)
tValue <- paste("t-value: ", round(test$statistic,4))
tpValue <- paste("p-value: ", round(test$p.value,4))
correlationValue <- paste("Correlation: ", round(test$estimate,4))
forest.perimeter.area.ratio=c(testMethod,testAlternative,tValue,
tpValue,correlationValue)
forest.perimeter.area.ratio
## [1] "Forest versus Perimeter/Area: Spearman's rank correlation rho"
## [2] "Alternative hypothesis: two.sided"
## [3] "t-value: 93056.4256"
## [4] "p-value: 0"
## [5] "Correlation: -0.9425"
############################################################
#
# save the results
#
############################################################
#write.csv(forest.metrics,file="forest-metrics.csv")
write(forest.metrics,file="forest-metrics.txt")
write(forest.core.correlation,file="forestcore-correlation.txt")
write(forest.perimeter,file="forest-perimeter.txt")
write(forest.perimeter.area.ratio,file="forest-ratio.txt")
############################################################
#
# Read Farm data
#
############################################################
farm = as.matrix(read.table("farm_1x1.txt"))
farm.r=raster(farm)
#Plotting farm
image(farm.r, main="farm", col=gray(1:0))
#Calculating the window size (Size of graph)
farm.lengthOfArea = ncol(farm)
farm.widthOfArea = nrow(farm)
############################################################
# farm output
############################################################
farm.total.area <- lsm_c_ca(farm.r, directions = 8)$value[2]*farm.lengthOfArea*farm.widthOfArea*10000
farm.prop.area <- lsm_c_ca(farm.r, directions = 8)$value[2]*10000
farm.number.patches <- lsm_c_np(farm.r, directions = 8)$value[2]
farm.total.edge <- lsm_l_te(farm.r, count_boundary = FALSE)$value/(10000)*farm.lengthOfArea*farm.widthOfArea*31.0461
farm.total.core.area <- lsm_c_tca(farm.r, directions = 8, consider_boundary = FALSE, edge_depth = 1)$value[2]*farm.lengthOfArea*farm.widthOfArea*10000
farm.prop.core.area <- lsm_c_tca(farm.r, directions = 8, consider_boundary = FALSE, edge_depth = 1)$value[2]*10000
farm.intermediate <- lsm_p_area(farm.r, directions = 8)
farm.max.patch.area <- max(farm.intermediate[which(farm.intermediate$class==1),]$value)*
farm.lengthOfArea*farm.widthOfArea*10000
farm.max.patch.index <- max(farm.intermediate[which(farm.intermediate$class==1),]$value)*10000
farm.mean.patch.area <- farm.total.area/farm.number.patches
farm.totalArea <- paste("Farm total area: ", farm.total.area)
farm.propTotalArea <- paste("Farm proportion total area: ", round(farm.prop.area,4))
farm.numberPatches <- paste("Farm number of patches: ", farm.number.patches)
farm.totalEdge <- paste("Farm total edge: ", round(farm.total.edge))
farm.totalCoreArea <- paste("Farm total core area: ", farm.total.core.area)
farm.propCoreArea <-paste("Farm proportion core area: ", round(farm.prop.core.area,4))
farm.maxPatchArea <- paste("Farm max patch area: ", farm.max.patch.area)
farm.maxPatchIndex <- paste("Farm max patch index: ", round(farm.max.patch.index,4))
farm.meanPatchArea <- paste("Farm mean patch area: ", round(farm.mean.patch.area,4))
farm.metrics=c(farm.totalArea,farm.propTotalArea,farm.numberPatches,
farm.totalEdge,farm.totalCoreArea,farm.propCoreArea,
farm.maxPatchArea,farm.maxPatchIndex,farm.meanPatchArea)
farm.metrics
## [1] "Farm total area: 13700" "Farm proportion total area: 0.1353"
## [3] "Farm number of patches: 25" "Farm total edge: 3787"
## [5] "Farm total core area: 10309" "Farm proportion core area: 0.1018"
## [7] "Farm max patch area: 5100" "Farm max patch index: 0.0504"
## [9] "Farm mean patch area: 548"
############################################################
# For farm versus core area
############################################################
f.intermediate <- lsm_p_area(farm.r, directions = 8)
farm.patch.area <- f.intermediate[which(f.intermediate$class==1),]$value*farm.lengthOfArea*farm.widthOfArea*10000
c.intermediate <- lsm_p_core(farm.r, directions = 8, consider_boundary = FALSE, edge_depth = 1)
farm.path.core.area <- c.intermediate[which(c.intermediate$class==1),]$value*farm.lengthOfArea*farm.widthOfArea*10000
plot(farm.patch.area,farm.path.core.area,main="farm: patch size area vs core area")
test<-cor.test(farm.patch.area,farm.path.core.area)
testMethod <- paste("Farm versus Core Area:", test$method)
testAlternative <- paste("Alternative hypothesis: ", test$alternative)
tDf <- paste("Degrees of freedom: ", test$parameter)
tValue <- paste("t-value: ", round(test$statistic,4))
tpValue <- paste("p-value: ", round(test$p.value,4))
correlationValue <- paste("Correlation: ", round(test$estimate,4))
farm.core.correlation=c(testMethod,testAlternative,tDf,
tValue,tpValue,correlationValue)
farm.core.correlation
## [1] "Farm versus Core Area: Pearson's product-moment correlation"
## [2] "Alternative hypothesis: two.sided"
## [3] "Degrees of freedom: 23"
## [4] "t-value: 67.8503"
## [5] "p-value: 0"
## [6] "Correlation: 0.9975"
############################################################
# For farm versus perimeter
############################################################
p.intermediate <- lsm_p_perim(farm.r, directions = 8)
farm.patch.perim <- p.intermediate[which(p.intermediate$class==1),]$value*312.5
plot(farm.patch.area,farm.patch.perim,main="farm: patch size area vs perimeter")
test<-cor.test(farm.patch.area,farm.patch.perim)
testMethod <- paste("Farm versus Perimeter Correlation:", test$method)
testAlternative <- paste("Alternative hypothesis: ", test$alternative)
tDf <- paste("Degrees of freedom: ", test$parameter)
tValue <- paste("t-value: ", round(test$statistic,4))
tpValue <- paste("p-value: ", round(test$p.value,4))
correlationValue <- paste("Correlation: ", round(test$estimate,4))
farm.perimeter=c(testMethod,testAlternative,tDf,
tValue,tpValue,correlationValue)
farm.perimeter
## [1] "Farm versus Perimeter Correlation: Pearson's product-moment correlation"
## [2] "Alternative hypothesis: two.sided"
## [3] "Degrees of freedom: 23"
## [4] "t-value: 16.4787"
## [5] "p-value: 0"
## [6] "Correlation: 0.9602"
############################################################
# For farm versus perimeter area ratio
############################################################
intermediate <- lsm_p_para(farm.r, directions = 8)
farm.patch.perim.area.ratio <- intermediate[which(intermediate$class==1),]$value
plot(farm.patch.area,farm.patch.perim.area.ratio,main="farm: patch size area vs perimeter area ratio")
test<-cor.test(farm.patch.area,farm.patch.perim.area.ratio,method="spearman")
## Warning in cor.test.default(farm.patch.area, farm.patch.perim.area.ratio, :
## Cannot compute exact p-value with ties
testMethod <- paste("Farm versus Perimeter/Area: ", test$method)
testAlternative <- paste("Alternative hypothesis: ", test$alternative)
tValue <- paste("t-value: ", round(test$statistic,4))
tpValue <- paste("p-value: ", round(test$p.value,4))
correlationValue <- paste("Correlation: ", round(test$estimate,4))
farm.perimeter.area.ratio=c(testMethod,testAlternative,tValue,
tpValue,correlationValue)
farm.perimeter.area.ratio
## [1] "Farm versus Perimeter/Area: Spearman's rank correlation rho"
## [2] "Alternative hypothesis: two.sided"
## [3] "t-value: 5062.7252"
## [4] "p-value: 0"
## [5] "Correlation: -0.9472"
############################################################
#
# save the results
#
############################################################
#write.csv(farm.metrics,file="farm-metrics.csv")
write(farm.metrics,file="farm-metrics.txt")
write(farm.core.correlation,file="farmcore-correlation.txt")
write(farm.perimeter,file="farm-perimeter.txt")
write(farm.perimeter.area.ratio,file="farm-ratio.txt")
############################################################
#
# Read Rural data
#
############################################################
rural=as.matrix(read.table("rural_1x1.txt"))
rural.r=raster(rural)
#Plotting rural
image(rural.r, main="rural", col=gray(1:0))
#Calculating the window size (Size of graph)
rural.lengthOfArea = ncol(rural)
rural.widthOfArea = nrow(rural)
############################################################
#rural output
############################################################
rural.total.area <- lsm_c_ca(rural.r, directions = 8)$value[2]*rural.lengthOfArea*rural.widthOfArea*10000
rural.prop.area <- lsm_c_ca(rural.r, directions = 8)$value[2]*10000
rural.number.patches <- lsm_c_np(rural.r, directions = 8)$value[2]
rural.total.edge <- lsm_l_te(rural.r, count_boundary = FALSE)$value/(10000)*rural.lengthOfArea*rural.widthOfArea*31.0461
rural.total.core.area <- lsm_c_tca(rural.r, directions = 8, consider_boundary = FALSE, edge_depth = 1)$value[2]*rural.lengthOfArea*rural.widthOfArea*10000
rural.prop.core.area <- lsm_c_tca(rural.r, directions = 8, consider_boundary = FALSE, edge_depth = 1)$value[2]*10000
rural.intermediate <- lsm_p_area(rural.r, directions = 8)
rural.max.patch.area <- max(rural.intermediate[which(rural.intermediate$class==1),]$value)*
rural.lengthOfArea*rural.widthOfArea*10000
rural.max.patch.index <- max(rural.intermediate[which(rural.intermediate$class==1),]$value)*10000
rural.mean.patch.area <- rural.total.area/rural.number.patches
rural.totalArea <- paste("Rural total area: ", rural.total.area)
rural.propTotalArea <- paste("Rural proportion total area: ", round(rural.prop.area,4))
rural.numberPatches <- paste("Rural number of patches: ", rural.number.patches)
rural.totalEdge <- paste("Rural total edge: ", round(rural.total.edge))
rural.totalCoreArea <- paste("Rural total core area: ", rural.total.core.area)
rural.propCoreArea <-paste("Rural proportion core area: ", round(rural.prop.core.area,4))
rural.maxPatchArea <- paste("Rural max patch area: ", rural.max.patch.area)
rural.maxPatchIndex <- paste("Rural max patch index: ", round(rural.max.patch.index,4))
rural.meanPatchArea <- paste("Rural mean patch area: ", round(rural.mean.patch.area,4))
rural.metrics=c(rural.totalArea,rural.propTotalArea,rural.numberPatches,
rural.totalEdge,rural.totalCoreArea,rural.propCoreArea,
rural.maxPatchArea,rural.maxPatchIndex,rural.meanPatchArea)
rural.metrics
## [1] "Rural total area: 12850"
## [2] "Rural proportion total area: 0.1269"
## [3] "Rural number of patches: 43"
## [4] "Rural total edge: 4680"
## [5] "Rural total core area: 8745"
## [6] "Rural proportion core area: 0.0864"
## [7] "Rural max patch area: 4675"
## [8] "Rural max patch index: 0.0462"
## [9] "Rural mean patch area: 298.8372"
############################################################
# For rural versus core area
############################################################
f.intermediate <- lsm_p_area(rural.r, directions = 8)
rural.patch.area <- f.intermediate[which(f.intermediate$class==1),]$value*rural.lengthOfArea*rural.widthOfArea*10000
c.intermediate <- lsm_p_core(rural.r, directions = 8, consider_boundary = FALSE, edge_depth = 1)
rural.path.core.area <- c.intermediate[which(c.intermediate$class==1),]$value*rural.lengthOfArea*rural.widthOfArea*10000
plot(rural.patch.area,rural.path.core.area,main="rural: patch size area vs core area")
test<-cor.test(rural.patch.area,rural.path.core.area)
testMethod <- paste("Rural versus Core Area:", test$method)
testAlternative <- paste("Alternative hypothesis: ", test$alternative)
tDf <- paste("Degrees of freedom: ", test$parameter)
tValue <- paste("t-value: ", round(test$statistic,4))
tpValue <- paste("p-value: ", round(test$p.value,4))
correlationValue <- paste("Correlation: ", round(test$estimate,4))
rural.core.correlation=c(testMethod,testAlternative,tDf,
tValue,tpValue,correlationValue)
rural.core.correlation
## [1] "Rural versus Core Area: Pearson's product-moment correlation"
## [2] "Alternative hypothesis: two.sided"
## [3] "Degrees of freedom: 41"
## [4] "t-value: 86.9125"
## [5] "p-value: 0"
## [6] "Correlation: 0.9973"
############################################################
# For rural versus perimeter
############################################################
p.intermediate <- lsm_p_perim(rural.r, directions = 8)
rural.patch.perim <- p.intermediate[which(p.intermediate$class==1),]$value*312.5
plot(rural.patch.area,rural.patch.perim,main="rural: patch size area vs perimeter")
test<-cor.test(rural.patch.area,rural.patch.perim)
testMethod <- paste("Rural versus Perimeter Correlation:", test$method)
testAlternative <- paste("Alternative hypothesis: ", test$alternative)
tDf <- paste("Degrees of freedom: ", test$parameter)
tValue <- paste("t-value: ", round(test$statistic,4))
tpValue <- paste("p-value: ", round(test$p.value,4))
correlationValue <- paste("Correlation: ", round(test$estimate,4))
rural.perimeter=c(testMethod,testAlternative,tDf,
tValue,tpValue,correlationValue)
rural.perimeter
## [1] "Rural versus Perimeter Correlation: Pearson's product-moment correlation"
## [2] "Alternative hypothesis: two.sided"
## [3] "Degrees of freedom: 41"
## [4] "t-value: 25.2179"
## [5] "p-value: 0"
## [6] "Correlation: 0.9692"
############################################################
# For rural versus perimeter area ratio
############################################################
intermediate <- lsm_p_para(rural.r, directions = 8)
rural.patch.perim.area.ratio <- intermediate[which(intermediate$class==1),]$value
plot(rural.patch.area,rural.patch.perim.area.ratio,main="rural: patch size area vs perimeter area ratio")
test<-cor.test(rural.patch.area,rural.patch.perim.area.ratio,method="spearman")
## Warning in cor.test.default(rural.patch.area, rural.patch.perim.area.ratio, :
## Cannot compute exact p-value with ties
testMethod <- paste("Rural versus Perimeter/Area: ", test$method)
testAlternative <- paste("Alternative hypothesis: ", test$alternative)
tValue <- paste("t-value: ", round(test$statistic,4))
tpValue <- paste("p-value: ", round(test$p.value,4))
correlationValue <- paste("Correlation: ", round(test$estimate,4))
rural.perimeter.area.ratio=c(testMethod,testAlternative,tValue,
tpValue,correlationValue)
rural.perimeter.area.ratio
## [1] "Rural versus Perimeter/Area: Spearman's rank correlation rho"
## [2] "Alternative hypothesis: two.sided"
## [3] "t-value: 25323.1199"
## [4] "p-value: 0"
## [5] "Correlation: -0.912"
############################################################
#
# save the results
#
############################################################
write(rural.metrics,file="rural-metrics.txt")
write(rural.core.correlation,file="ruralcore-correlation.txt")
write(rural.perimeter,file="rural-perimeter.txt")
write(rural.perimeter.area.ratio,file="rural-ratio.txt")
############################################################
#
# Read Urban data
#
############################################################
urban=as.matrix(read.table("urban_1x1.txt"))
urban.r=raster(urban)
#Plotting urban
image(urban.r, main="urban", col=gray(1:0))
#Calculating the window size (Size of graph)
urban.lengthOfArea = ncol(urban)
urban.widthOfArea = nrow(urban)
############################################################
# urban output
############################################################
urban.total.area <- lsm_c_ca(urban.r, directions = 8)$value[2]*urban.lengthOfArea*urban.widthOfArea*10000
urban.prop.area <- lsm_c_ca(urban.r, directions = 8)$value[2]*10000
urban.number.patches <- lsm_c_np(urban.r, directions = 8)$value[2]
urban.total.edge <- lsm_l_te(urban.r, count_boundary = FALSE)$value/(10000)*urban.lengthOfArea*urban.widthOfArea*31.0461
urban.total.core.area <- lsm_c_tca(urban.r, directions = 8, consider_boundary = FALSE, edge_depth = 1)$value[2]*urban.lengthOfArea*urban.widthOfArea*10000
urban.prop.core.area <- lsm_c_tca(urban.r, directions = 8, consider_boundary = FALSE, edge_depth = 1)$value[2]*10000
urban.intermediate <- lsm_p_area(urban.r, directions = 8)
urban.max.patch.area <- max(urban.intermediate[which(urban.intermediate$class==1),]$value)*
urban.lengthOfArea*urban.widthOfArea*10000
urban.max.patch.index <- max(urban.intermediate[which(urban.intermediate$class==1),]$value)*10000
urban.mean.patch.area <- urban.total.area/urban.number.patches
urban.totalArea <- paste("Urban total area: ", urban.total.area)
urban.propTotalArea <- paste("Urban proportion total area: ", round(urban.prop.area,4))
urban.numberPatches <- paste("Urban number of patches: ", urban.number.patches)
urban.totalEdge <- paste("Urban total edge: ", round(urban.total.edge))
urban.totalCoreArea <- paste("Urban total core area: ", urban.total.core.area)
urban.propCoreArea <-paste("Urban proportion core area: ", round(urban.prop.core.area,4))
urban.maxPatchArea <- paste("Urban max patch area: ", urban.max.patch.area)
urban.maxPatchIndex <- paste("Urban max patch index: ", round(urban.max.patch.index,4))
urban.meanPatchArea <- paste("Urban mean patch area: ", round(urban.mean.patch.area,4))
urban.metrics=c(urban.totalArea,urban.propTotalArea,urban.numberPatches,
urban.totalEdge,urban.totalCoreArea,urban.propCoreArea,
urban.maxPatchArea,urban.maxPatchIndex,urban.meanPatchArea)
urban.metrics
## [1] "Urban total area: 21675"
## [2] "Urban proportion total area: 0.2141"
## [3] "Urban number of patches: 28"
## [4] "Urban total edge: 4833"
## [5] "Urban total core area: 17292"
## [6] "Urban proportion core area: 0.1708"
## [7] "Urban max patch area: 13225"
## [8] "Urban max patch index: 0.1306"
## [9] "Urban mean patch area: 774.1071"
############################################################
# For urban versus core area
############################################################
f.intermediate <- lsm_p_area(urban.r, directions = 8)
urban.patch.area <- f.intermediate[which(f.intermediate$class==1),]$value*urban.lengthOfArea*urban.widthOfArea*10000
c.intermediate <- lsm_p_core(urban.r, directions = 8, consider_boundary = FALSE, edge_depth = 1)
urban.path.core.area <- c.intermediate[which(c.intermediate$class==1),]$value*urban.lengthOfArea*urban.widthOfArea*10000
plot(urban.patch.area,urban.path.core.area,main="urban: patch size area vs core area")
test<-cor.test(urban.patch.area,urban.path.core.area)
testMethod <- paste("Urban versus Core Area:", test$method)
testAlternative <- paste("Alternative hypothesis: ", test$alternative)
tDf <- paste("Degrees of freedom: ", test$parameter)
tValue <- paste("t-value: ", round(test$statistic,4))
tpValue <- paste("p-value: ", round(test$p.value,4))
correlationValue <- paste("Correlation: ", round(test$estimate,4))
urban.core.correlation=c(testMethod,testAlternative,tDf,
tValue,tpValue,correlationValue)
urban.core.correlation
## [1] "Urban versus Core Area: Pearson's product-moment correlation"
## [2] "Alternative hypothesis: two.sided"
## [3] "Degrees of freedom: 26"
## [4] "t-value: 854.1628"
## [5] "p-value: 0"
## [6] "Correlation: 1"
############################################################
#For urban versus perimeter
############################################################
p.intermediate <- lsm_p_perim(urban.r, directions = 8)
urban.patch.perim <- p.intermediate[which(p.intermediate$class==1),]$value*312.5
plot(urban.patch.area,urban.patch.perim,main="urban: patch size area vs perimeter")
test<-cor.test(urban.patch.area,urban.patch.perim)
testMethod <- paste("Urban versus Perimeter Correlation:", test$method)
testAlternative <- paste("Alternative hypothesis: ", test$alternative)
tDf <- paste("Degrees of freedom: ", test$parameter)
tValue <- paste("t-value: ", round(test$statistic,4))
tpValue <- paste("p-value: ", round(test$p.value,4))
correlationValue <- paste("Correlation: ", round(test$estimate,4))
urban.perimeter=c(testMethod,testAlternative,tDf,
tValue,tpValue,correlationValue)
urban.perimeter
## [1] "Urban versus Perimeter Correlation: Pearson's product-moment correlation"
## [2] "Alternative hypothesis: two.sided"
## [3] "Degrees of freedom: 26"
## [4] "t-value: 238.4164"
## [5] "p-value: 0"
## [6] "Correlation: 0.9998"
############################################################
#For urban versus perimeter area ratio
############################################################
intermediate <- lsm_p_para(urban.r, directions = 8)
urban.patch.perim.area.ratio <- intermediate[which(intermediate$class==1),]$value
plot(urban.patch.area,urban.patch.perim.area.ratio,main="urban: patch size area vs perimeter area ratio")
test<-cor.test(urban.patch.area,urban.patch.perim.area.ratio,method="spearman")
## Warning in cor.test.default(urban.patch.area, urban.patch.perim.area.ratio, :
## Cannot compute exact p-value with ties
testMethod <- paste("Urban versus Perimeter/Area: ", test$method)
testAlternative <- paste("Alternative hypothesis: ", test$alternative)
tValue <- paste("t-value: ", round(test$statistic,4))
tpValue <- paste("p-value: ", round(test$p.value,4))
correlationValue <- paste("Correlation: ", round(test$estimate,4))
urban.perimeter.area.ratio=c(testMethod,testAlternative,tValue,
tpValue,correlationValue)
urban.perimeter.area.ratio
## [1] "Urban versus Perimeter/Area: Spearman's rank correlation rho"
## [2] "Alternative hypothesis: two.sided"
## [3] "t-value: 7135.1845"
## [4] "p-value: 0"
## [5] "Correlation: -0.9527"
############################################################
#
# save the results
#
############################################################
write(urban.metrics,file="urban-metrics.txt")
write(urban.core.correlation,file="urbancore-correlation.txt")
write(urban.perimeter,file="urban-perimeter.txt")
write(urban.perimeter.area.ratio,file="urban-ratio.txt")
############################################################