suppressWarnings({library(sp)
library(readxl)
library(dplyr)
library(glmnet)
library(car)
library(broom)
library(estimatr)
library(lm.beta)
library(tidyverse)
})## The legacy packages maptools, rgdal, and rgeos, underpinning the sp package,
## which was just loaded, will retire in October 2023.
## Please refer to R-spatial evolution reports for details, especially
## https://r-spatial.org/r/2023/05/15/evolution4.html.
## It may be desirable to make the sf package available;
## package maintainers should consider adding sf to Suggests:.
## The sp package is now running under evolution status 2
## (status 2 uses the sf package in place of rgdal)##
## Attaching package: 'dplyr'## The following objects are masked from 'package:stats':
##
## filter, lag## The following objects are masked from 'package:base':
##
## intersect, setdiff, setequal, union## Loading required package: Matrix## Loaded glmnet 4.1-7## Loading required package: carData##
## Attaching package: 'car'## The following object is masked from 'package:dplyr':
##
## recode## ── Attaching core tidyverse packages ──────────────────────── tidyverse 2.0.0 ──
## ✔ forcats 1.0.0 ✔ readr 2.1.4
## ✔ ggplot2 3.4.2 ✔ stringr 1.5.0
## ✔ lubridate 1.9.2 ✔ tibble 3.2.1
## ✔ purrr 1.0.1 ✔ tidyr 1.3.0
## ── Conflicts ────────────────────────────────────────── tidyverse_conflicts() ──
## ✖ tidyr::expand() masks Matrix::expand()
## ✖ dplyr::filter() masks stats::filter()
## ✖ dplyr::lag() masks stats::lag()
## ✖ tidyr::pack() masks Matrix::pack()
## ✖ car::recode() masks dplyr::recode()
## ✖ purrr::some() masks car::some()
## ✖ tidyr::unpack() masks Matrix::unpack()
## ℹ Use the conflicted package (<http://conflicted.r-lib.org/>) to force all conflicts to become errorsAMAZON <- read_xlsx("M:/lvillaquiran3/CAPSTONE/PRODES/CAPSTONE_LRM/RMS_POINTS_Merg3_TableToExcel.xlsx")
head(AMAZON)## # A tibble: 6 × 6
## FID CID Distance_M Distance_I Distance_R Slope__Deg
## <dbl> <dbl> <dbl> <dbl> <dbl> <dbl>
## 1 0 0 328164 236946 189942 0
## 2 1 0 522987 140502 183293 0.0573
## 3 2 0 155108 57099. 213652 0
## 4 3 0 109235 29802. 77912 0
## 5 4 0 309279 265613 271354 0
## 6 5 0 43322. 17693. 48991. 0dependent.variable <- AMAZON$CID
independent.var <- AMAZON %>% select(Distance_M, Distance_I, Distance_R, Slope__Deg)
independent.var3 <- AMAZON %>% select(Distance_M, Distance_I, Distance_R)
independent.var2 <- AMAZON %>% select(Distance_R, Distance_I)
independent.var1 <- AMAZON %>% select(Distance_R)OGmodel <- glm(dependent.variable ~., data = cbind(dependent.variable, independent.var), family = "binomial") ## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurredmodel.sum <- summary(OGmodel)
print(model.sum)##
## Call:
## glm(formula = dependent.variable ~ ., family = "binomial", data = cbind(dependent.variable,
## independent.var))
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) 9.610e-01 3.354e-01 2.866 0.004163 **
## Distance_M -1.431e-06 1.464e-06 -0.977 0.328372
## Distance_I 2.063e-05 6.030e-06 3.421 0.000623 ***
## Distance_R -2.622e-05 6.204e-06 -4.226 2.38e-05 ***
## Slope__Deg 2.671e-01 3.836e-01 0.696 0.486289
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for binomial family taken to be 1)
##
## Null deviance: 277.26 on 199 degrees of freedom
## Residual deviance: 245.43 on 195 degrees of freedom
## AIC: 255.43
##
## Number of Fisher Scoring iterations: 12anova1<- anova(OGmodel, test = 'Chisq')
print(anova1)## Analysis of Deviance Table
##
## Model: binomial, link: logit
##
## Response: dependent.variable
##
## Terms added sequentially (first to last)
##
##
## Df Deviance Resid. Df Resid. Dev Pr(>Chi)
## NULL 199 277.26
## Distance_M 1 2.6170 198 274.64 0.10573
## Distance_I 1 3.2715 197 271.37 0.07049 .
## Distance_R 1 24.6103 196 246.76 7.018e-07 ***
## Slope__Deg 1 1.3292 195 245.43 0.24894
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1model3 <- glm(dependent.variable ~., data = cbind(dependent.variable, independent.var3), family = "binomial")
model.sum3 <- summary(model3)
print(model.sum3)##
## Call:
## glm(formula = dependent.variable ~ ., family = "binomial", data = cbind(dependent.variable,
## independent.var3))
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) 1.064e+00 3.173e-01 3.353 0.000800 ***
## Distance_M -1.685e-06 1.442e-06 -1.168 0.242644
## Distance_I 2.102e-05 6.021e-06 3.491 0.000482 ***
## Distance_R -2.671e-05 6.190e-06 -4.315 1.6e-05 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for binomial family taken to be 1)
##
## Null deviance: 277.26 on 199 degrees of freedom
## Residual deviance: 246.76 on 196 degrees of freedom
## AIC: 254.76
##
## Number of Fisher Scoring iterations: 4model2 <- glm(dependent.variable ~., data = cbind(dependent.variable, independent.var2), family = "binomial")
model.sum2 <- summary(model2)
print(model.sum2)##
## Call:
## glm(formula = dependent.variable ~ ., family = "binomial", data = cbind(dependent.variable,
## independent.var2))
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) 8.121e-01 2.296e-01 3.537 0.000405 ***
## Distance_R -2.538e-05 6.038e-06 -4.203 2.63e-05 ***
## Distance_I 1.839e-05 5.536e-06 3.322 0.000894 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for binomial family taken to be 1)
##
## Null deviance: 277.26 on 199 degrees of freedom
## Residual deviance: 248.15 on 197 degrees of freedom
## AIC: 254.15
##
## Number of Fisher Scoring iterations: 4model1 <- glm(dependent.variable ~., data = cbind(dependent.variable, independent.var1), family = "binomial")
model.sum1 <- summary(model1)
print(model.sum1)##
## Call:
## glm(formula = dependent.variable ~ ., family = "binomial", data = cbind(dependent.variable,
## independent.var1))
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) 6.107e-01 2.143e-01 2.850 0.004376 **
## Distance_R -6.860e-06 1.841e-06 -3.726 0.000194 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for binomial family taken to be 1)
##
## Null deviance: 277.26 on 199 degrees of freedom
## Residual deviance: 261.45 on 198 degrees of freedom
## AIC: 265.45
##
## Number of Fisher Scoring iterations: 4anova2<- anova(model2, test = 'Chisq')
print(anova2)## Analysis of Deviance Table
##
## Model: binomial, link: logit
##
## Response: dependent.variable
##
## Terms added sequentially (first to last)
##
##
## Df Deviance Resid. Df Resid. Dev Pr(>Chi)
## NULL 199 277.26
## Distance_R 1 15.804 198 261.45 7.025e-05 ***
## Distance_I 1 13.301 197 248.15 0.0002653 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1AIC_TABLE <- data.frame(
"Model" = c("M1", "M2", "M3", "M4"),
"Variables" = c("Illegal Mines, Indigenous Territories, Roads, slope", "Illegal Mines, Indigenous Territories, Roads", "Indigenous Territory, Roads", "Roads"),
"AIC" = c(255.43
,254.76, 254.15, 265.45)
)
# Print the data frame
print(AIC_TABLE)## Model Variables AIC
## 1 M1 Illegal Mines, Indigenous Territories, Roads, slope 255.43
## 2 M2 Illegal Mines, Indigenous Territories, Roads 254.76
## 3 M3 Indigenous Territory, Roads 254.15
## 4 M4 Roads 265.45library(flextable)## Warning: package 'flextable' was built under R version 4.3.1##
## Attaching package: 'flextable'## The following object is masked from 'package:purrr':
##
## composeflextable(AIC_TABLE)Model | Variables | AIC |
|---|---|---|
M1 | Illegal Mines, Indigenous Territories, Roads, slope | 255.43 |
M2 | Illegal Mines, Indigenous Territories, Roads | 254.76 |
M3 | Indigenous Territory, Roads | 254.15 |
M4 | Roads | 265.45 |
BETA_TABLE <- data.frame(
"Variable" = c("Distance_M", "Distance_I", "Distance_R", "Slope_Des"),
"Beta_coef" = c(-1.43101534945699e-06, 2.06310332743348e-05,-2.62177148076232e-05, 0.267069903487502))
print(BETA_TABLE)## Variable Beta_coef
## 1 Distance_M -1.431015e-06
## 2 Distance_I 2.063103e-05
## 3 Distance_R -2.621771e-05
## 4 Slope_Des 2.670699e-01write.csv(BETA_TABLE,file = "BETA_TABLE.txt", row.names = TRUE)
flextable(BETA_TABLE)Variable | Beta_coef |
|---|---|
Distance_M | -0.000001431015 |
Distance_I | 0.000020631033 |
Distance_R | -0.000026217715 |
Slope_Des | 0.267069903488 |