Loading the libraries to be used
library(readxl)
library(dplyr)
library(lme4)
library(sjPlot)
library(psycho)
library(tidyverse)
library(lmerTest)rm(list=ls(all=TRUE))
getwd()
setwd ("C:/Users/moliv/Documents/R")#reading data sets from csv files
btex1<-read.csv("BTEX_1.csv", header=T, sep=",")
btex2<-read.csv("BTEX_2.csv", header=T, sep=",")
btex3<-read.csv("BTEX_3.csv", header=T, sep=",")vasilikos_stations <- c("Vasilikos Area- Gov. Beach","Vasilikos Area- Zygi")btex1[,c(2:6)] <- sapply(btex1[,c(2:6)], as.numeric)
btex1$WEEK <- as.factor("1")
btex1$ID <- as.factor(as.integer(btex1$ID))
btex1$LOCATION <- ifelse(btex1$Address %in% vasilikos_stations,"1","0")
btex1$LOCATION <- as.factor(btex1$LOCATION)
str(btex1)btex2[,c(2:6)] <- sapply(btex1[,c(2:6)], as.numeric)## Warning in matrix(value, n, p): data length [215] is not a sub-multiple or
## multiple of the number of rows [41]btex2$WEEK <- as.factor("2")
btex2$ID <- as.factor(as.integer(btex2$ID))
btex2$LOCATION <- ifelse(btex2$Address %in% vasilikos_stations,"1","0")
btex2$LOCATION <- as.factor(btex2$LOCATION)
str(btex2)btex3[,c(2:6)] <- sapply(btex1[,c(2:6)], as.numeric)
btex3$WEEK <- as.factor("3")
btex3$ID <- as.factor(as.integer(btex3$ID))
btex3$LOCATION <- ifelse(btex3$Address %in% vasilikos_stations,"1","0")
btex3$LOCATION <- as.factor(btex3$LOCATION)
str(btex3)btex <- rbind(btex1,btex2,btex3)
str(btex)## 'data.frame': 127 obs. of 11 variables:
## $ CODE_ID : Factor w/ 127 levels "1_1","10_1","11_1",..: 1 12 22 32 39 40 41 42 43 2 ...
## $ B78 : num 0.0298 0.0895 0.0895 0.1836 0.0298 ...
## $ T91 : num 0.159 0.159 0.529 0.473 0.159 ...
## $ EB91 : num 0.0161 0.0161 0.0161 0.0161 0.0161 ...
## $ PMX106 : num 0.0345 0.1035 0.3943 0.2362 0.1035 ...
## $ OX106 : num 0.24 0.368 1.029 0.909 0.309 ...
## $ WEEK : Factor w/ 3 levels "1","2","3": 1 1 1 1 1 1 1 1 1 1 ...
## $ ID : Factor w/ 45 levels "1","2","3","4",..: 1 2 3 4 5 6 7 8 9 10 ...
## $ GPS : Factor w/ 45 levels "34°40'46.2\"N 33°02'38.3\"E",..: 32 38 35 36 40 42 43 41 26 29 ...
## $ Address : Factor w/ 13 levels "Asgata","Choirokoitia, Cyprus",..: 10 10 10 10 2 2 2 2 9 9 ...
## $ LOCATION: Factor w/ 2 levels "0","1": 1 1 1 1 1 1 1 1 1 1 ...knitr::kable(head(btex))| CODE_ID | B78 | T91 | EB91 | PMX106 | OX106 | WEEK | ID | GPS | Address | LOCATION |
|---|---|---|---|---|---|---|---|---|---|---|
| 1_1 | 0.0298498 | 0.1592351 | 0.0160969 | 0.0345077 | 0.2396894 | 1 | 1 | 34°46’33.0“N 33°19’19.9”E | Tochni, Cyprus | 0 |
| 2_1 | 0.0895493 | 0.1592351 | 0.0160969 | 0.1035230 | 0.3684043 | 1 | 2 | 34°47’02.6“N 33°19’18.0”E | Tochni, Cyprus | 0 |
| 3_1 | 0.0895493 | 0.5286092 | 0.0160969 | 0.3942761 | 1.0285316 | 1 | 3 | 34°46’48.2“N 33°19’11.6”E | Tochni, Cyprus | 0 |
| 4_1 | 0.1836340 | 0.4734591 | 0.0160969 | 0.2362432 | 0.9094332 | 1 | 4 | 34°46’48.8“N 33°19’33.4”E | Tochni, Cyprus | 0 |
| 5_1 | 0.0298498 | 0.1592351 | 0.0160969 | 0.1035230 | 0.3087047 | 1 | 5 | 34°47’39.4“N 33°20’12.1”E | Choirokoitia, Cyprus | 0 |
| 6_1 | 0.0895493 | 0.4582835 | 0.0160969 | 0.2195806 | 0.7835104 | 1 | 6 | 34°47’56.1“N 33°20’00.6”E | Choirokoitia, Cyprus | 0 |
print(paste("We have",length(levels(btex$Address)),"different Communities including the two Vasilikos stations"))## [1] "We have 13 different Communities including the two Vasilikos stations"knitr::kable(levels(btex$Address), caption = "Communities")| x |
|---|
| Asgata |
| Choirokoitia, Cyprus |
| Kalavasos |
| Limassol |
| Limassol Govt Air Station |
| Mari, Cyprus |
| Maroni, Cyprus |
| Pentakomo |
| Psematismenos, Cyprus |
| Tochni, Cyprus |
| Vasilikos Area- Gov. Beach |
| Vasilikos Area- Zygi |
| Zygi |
print(paste("We have",length(levels(btex$GPS)),"different sample points"))## [1] "We have 45 different sample points"Using GPS Coordinates (Sample Points) as Random Effect Location = 1 for Vasilikos Stations (“Vasilikos Area- Gov. Beach” and “Vasilikos Area- Zygi”) and 0 for No Vasilikos stations
#model <- lmer(B78~week + location + (1|GPS), data=btex)
model2 <- lmer(B78~WEEK + LOCATION + Address + (1|GPS), data=btex)## fixed-effect model matrix is rank deficient so dropping 1 column / coefficient## fixed-effect model matrix is rank deficient so dropping 1 column / coefficient
model3 <- lmer(T91~WEEK + LOCATION + Address + (1|GPS), data=btex)## fixed-effect model matrix is rank deficient so dropping 1 column / coefficient## boundary (singular) fit: see ?isSingular## fixed-effect model matrix is rank deficient so dropping 1 column / coefficient
model4 <- lmer(EB91~WEEK + LOCATION + Address + (1|GPS), data=btex)## fixed-effect model matrix is rank deficient so dropping 1 column / coefficient## boundary (singular) fit: see ?isSingular## fixed-effect model matrix is rank deficient so dropping 1 column / coefficient
model5 <- lmer(PMX106~WEEK + LOCATION + Address + (1|GPS), data=btex)## fixed-effect model matrix is rank deficient so dropping 1 column / coefficient## boundary (singular) fit: see ?isSingular## fixed-effect model matrix is rank deficient so dropping 1 column / coefficient
model6 <- lmer(OX106~WEEK + LOCATION + Address + (1|GPS), data=btex)## fixed-effect model matrix is rank deficient so dropping 1 column / coefficient## boundary (singular) fit: see ?isSingular## fixed-effect model matrix is rank deficient so dropping 1 column / coefficient
models <- c(model2, model3, model4, model5, model6)
tab_model(models)| B78 | T91 | EB91 | PMX106 | OX106 | |||||||||||
|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|
| Predictors | Estimates | CI | p | Estimates | CI | p | Estimates | CI | p | Estimates | CI | p | Estimates | CI | p |
| (Intercept) | 0.43 | 0.21 – 0.65 | <0.001 | 2.55 | 1.72 – 3.39 | <0.001 | 0.35 | 0.14 – 0.56 | 0.001 | 0.65 | 0.06 – 1.24 | 0.031 | 1.88 | 0.09 – 3.67 | 0.040 |
| WEEK [2] | -0.01 | -0.13 – 0.12 | 0.933 | -0.06 | -0.61 – 0.49 | 0.833 | 0.07 | -0.07 – 0.21 | 0.345 | -0.12 | -0.51 – 0.27 | 0.544 | -0.40 | -1.58 – 0.78 | 0.509 |
| WEEK [3] | 0.00 | -0.12 – 0.13 | 0.944 | -0.02 | -0.56 – 0.53 | 0.954 | 0.00 | -0.14 – 0.14 | 0.982 | 0.03 | -0.36 – 0.41 | 0.899 | 0.09 | -1.08 – 1.26 | 0.878 |
| LOCATION [1] | -0.15 | -0.47 – 0.16 | 0.343 | -1.43 | -2.60 – -0.25 | 0.017 | -0.26 | -0.56 – 0.04 | 0.086 | 0.01 | -0.82 – 0.84 | 0.983 | -0.01 | -2.53 – 2.51 | 0.994 |
|
Address [Choirokoitia, Cyprus] |
-0.36 | -0.64 – -0.07 | 0.014 | -2.21 | -3.26 – -1.16 | <0.001 | -0.36 | -0.63 – -0.09 | 0.009 | -0.50 | -1.25 – 0.25 | 0.188 | -1.28 | -3.54 – 0.97 | 0.265 |
| Address [Kalavasos] | 0.10 | -0.19 – 0.38 | 0.498 | -0.65 | -1.70 – 0.40 | 0.227 | -0.01 | -0.28 – 0.26 | 0.931 | 1.54 | 0.79 – 2.28 | <0.001 | 3.76 | 1.50 – 6.01 | 0.001 |
| Address [Limassol] | 1.11 | 0.80 – 1.42 | <0.001 | 1.49 | 0.35 – 2.62 | 0.010 | 0.52 | 0.23 – 0.81 | <0.001 | 2.59 | 1.78 – 3.39 | <0.001 | 7.60 | 5.17 – 10.03 | <0.001 |
|
Address [Limassol Govt Air Station] |
0.44 | -0.07 – 0.95 | 0.091 | -0.02 | -1.96 – 1.93 | 0.987 | 0.16 | -0.34 – 0.65 | 0.533 | 1.47 | 0.09 – 2.85 | 0.037 | 5.15 | 0.98 – 9.32 | 0.016 |
| Address [Mari, Cyprus] | -0.26 | -0.55 – 0.03 | 0.074 | -1.78 | -2.84 – -0.73 | 0.001 | -0.31 | -0.58 – -0.04 | 0.022 | -0.29 | -1.03 – 0.46 | 0.453 | -0.77 | -3.03 – 1.49 | 0.503 |
| Address [Maroni, Cyprus] | -0.21 | -0.50 – 0.08 | 0.150 | -2.04 | -3.12 – -0.97 | <0.001 | -0.35 | -0.63 – -0.08 | 0.011 | -0.44 | -1.20 – 0.33 | 0.263 | -0.96 | -3.27 – 1.35 | 0.415 |
| Address [Pentakomo] | -0.24 | -0.53 – 0.05 | 0.107 | -1.19 | -2.27 – -0.12 | 0.030 | -0.26 | -0.54 – 0.01 | 0.062 | -0.28 | -1.04 – 0.48 | 0.467 | 1.28 | -1.02 – 3.59 | 0.275 |
|
Address [Psematismenos, Cyprus] |
-0.32 | -0.61 – -0.03 | 0.030 | -2.21 | -3.29 – -1.14 | <0.001 | -0.36 | -0.63 – -0.09 | 0.010 | -0.46 | -1.22 – 0.30 | 0.235 | -1.17 | -3.47 – 1.14 | 0.322 |
| Address [Tochni, Cyprus] | -0.33 | -0.62 – -0.05 | 0.022 | -2.24 | -3.30 – -1.19 | <0.001 | -0.08 | -0.35 – 0.18 | 0.537 | -0.48 | -1.22 – 0.27 | 0.209 | -1.01 | -3.27 – 1.24 | 0.379 |
|
Address [Vasilikos Area- Gov. Beach] |
0.01 | -0.36 – 0.38 | 0.954 | -0.06 | -1.42 – 1.31 | 0.935 | 0.01 | -0.34 – 0.35 | 0.975 | 0.16 | -0.80 – 1.13 | 0.738 | 0.62 | -2.30 – 3.55 | 0.675 |
| Address [Zygi] | -0.10 | -0.40 – 0.19 | 0.501 | -1.54 | -2.65 – -0.44 | 0.006 | -0.27 | -0.55 – 0.01 | 0.059 | -0.07 | -0.85 – 0.71 | 0.854 | 0.03 | -2.34 – 2.39 | 0.982 |
| Random Effects | |||||||||||||||
| σ2 | 0.09 | 1.66 | 0.11 | 0.83 | 7.61 | ||||||||||
| τ00 | 0.01 GPS | 0.00 GPS | 0.00 GPS | 0.00 GPS | 0.00 GPS | ||||||||||
| ICC | 0.11 | ||||||||||||||
| N | 45 GPS | 45 GPS | 45 GPS | 45 GPS | 45 GPS | ||||||||||
| Observations | 127 | 127 | 127 | 127 | 127 | ||||||||||
| Marginal R2 / Conditional R2 | 0.573 / 0.621 | 0.395 / NA | 0.344 / NA | 0.497 / NA | 0.448 / NA | ||||||||||