library(readxl)
locations <- read_excel("locations.xlsx")
To make maps on R, you need to register with Google Cloud Platform (free trial) and get an API key. Walkthrough for setting this up is here -> https://www.littlemissdata.com/blog/maps
Please note that the above walkthrough does not tell you to register Geocoding, only API Maps Static. You need both to do what I have done here. It will give you one key, with both registered to it.
You need to load ggmap and register your google API key for API Maps Static and Geocoding. You also need an internet connection to knit it as it is using google maps to find the coordinates.
library(ggplot2)
library(dplyr)
##
## 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
library(ggmap)
## Google's Terms of Service: https://cloud.google.com/maps-platform/terms/.
## Please cite ggmap if you use it! See citation("ggmap") for details.
ggmap::register_google(key = "AIzaSyB-MtKgp2mleNTDL2te3Y-JNBSX2klxVY0") #your API for Maps Static and Geocoding on Google Cloud. It must be set to unrestricted in Google.
geo <- geocode("Sydney", size = 0.3) #find coordinates of where you want the centre of your map to be.
## Source : https://maps.googleapis.com/maps/api/geocode/json?address=Sydney&key=xxx-MtKgp2mleNTDL2te3Y-JNBSX2klxVY0
print(geo)
## # A tibble: 1 x 2
## lon lat
## <dbl> <dbl>
## 1 151. -33.9
register_google(key = "AIzaSyB-MtKgp2mleNTDL2te3Y-JNBSX2klxVY0")
sc <- mutate_geocode(locations, SC) #calculate coordinates for each location. syntax = mutate_geocode(#YourDatset, #Your Locations Variable Name)
## Source : https://maps.googleapis.com/maps/api/geocode/json?address=Westfield+Eastgardens&key=xxx-MtKgp2mleNTDL2te3Y-JNBSX2klxVY0
## Source : https://maps.googleapis.com/maps/api/geocode/json?address=Fashion+Spree+Orange+Grove&key=xxx-MtKgp2mleNTDL2te3Y-JNBSX2klxVY0
## Source : https://maps.googleapis.com/maps/api/geocode/json?address=Neeta+City&key=xxx-MtKgp2mleNTDL2te3Y-JNBSX2klxVY0
## Source : https://maps.googleapis.com/maps/api/geocode/json?address=Marrickville+Metro&key=xxx-MtKgp2mleNTDL2te3Y-JNBSX2klxVY0
## Source : https://maps.googleapis.com/maps/api/geocode/json?address=Westfield+Penrith&key=xxx-MtKgp2mleNTDL2te3Y-JNBSX2klxVY0
## Source : https://maps.googleapis.com/maps/api/geocode/json?address=Westpoint+Shopping+Centre,+Blacktown&key=xxx-MtKgp2mleNTDL2te3Y-JNBSX2klxVY0
## Source : https://maps.googleapis.com/maps/api/geocode/json?address=Randwick+Shopping+Centre&key=xxx-MtKgp2mleNTDL2te3Y-JNBSX2klxVY0
## Source : https://maps.googleapis.com/maps/api/geocode/json?address=Kellyville+Village+Shops&key=xxx-MtKgp2mleNTDL2te3Y-JNBSX2klxVY0
## Source : https://maps.googleapis.com/maps/api/geocode/json?address=Rouse+Hill+Town+Centre&key=xxx-MtKgp2mleNTDL2te3Y-JNBSX2klxVY0
## Source : https://maps.googleapis.com/maps/api/geocode/json?address=Hurstville+Westfield&key=xxx-MtKgp2mleNTDL2te3Y-JNBSX2klxVY0
## Source : https://maps.googleapis.com/maps/api/geocode/json?address=Revesby+Village+Centre&key=xxx-MtKgp2mleNTDL2te3Y-JNBSX2klxVY0
## Source : https://maps.googleapis.com/maps/api/geocode/json?address=1-3+Brady+Street,+Mosman+NSW+2088&key=xxx-MtKgp2mleNTDL2te3Y-JNBSX2klxVY0
## Source : https://maps.googleapis.com/maps/api/geocode/json?address=Birkenhead+Point&key=xxx-MtKgp2mleNTDL2te3Y-JNBSX2klxVY0
## Source : https://maps.googleapis.com/maps/api/geocode/json?address=Westfield+Hornsby&key=xxx-MtKgp2mleNTDL2te3Y-JNBSX2klxVY0
## Source : https://maps.googleapis.com/maps/api/geocode/json?address=Minto+Marketplace&key=xxx-MtKgp2mleNTDL2te3Y-JNBSX2klxVY0
## Source : https://maps.googleapis.com/maps/api/geocode/json?address=Castle+Towers+Shopping+Centre,+6-14+Castle+Street,+Castle+Hill+NSW+2154&key=xxx-MtKgp2mleNTDL2te3Y-JNBSX2klxVY0
## Source : https://maps.googleapis.com/maps/api/geocode/json?address=Westfield+Parramatta&key=xxx-MtKgp2mleNTDL2te3Y-JNBSX2klxVY0
## Source : https://maps.googleapis.com/maps/api/geocode/json?address=Westfield+Bondi+Junction&key=xxx-MtKgp2mleNTDL2te3Y-JNBSX2klxVY0
## Source : https://maps.googleapis.com/maps/api/geocode/json?address=Top+Ryde+City+Shopping+Centre&key=xxx-MtKgp2mleNTDL2te3Y-JNBSX2klxVY0
## Source : https://maps.googleapis.com/maps/api/geocode/json?address=Rhodes+Waterside&key=xxx-MtKgp2mleNTDL2te3Y-JNBSX2klxVY0
## Source : https://maps.googleapis.com/maps/api/geocode/json?address=Westfield+Burwood&key=xxx-MtKgp2mleNTDL2te3Y-JNBSX2klxVY0
## Source : https://maps.googleapis.com/maps/api/geocode/json?address=DFO+Homebush&key=xxx-MtKgp2mleNTDL2te3Y-JNBSX2klxVY0
## Source : https://maps.googleapis.com/maps/api/geocode/json?address=Bankstown+Centro&key=xxx-MtKgp2mleNTDL2te3Y-JNBSX2klxVY0
## Source : https://maps.googleapis.com/maps/api/geocode/json?address=Ashfield+Mall&key=xxx-MtKgp2mleNTDL2te3Y-JNBSX2klxVY0
## Source : https://maps.googleapis.com/maps/api/geocode/json?address=Macquarie+Centre&key=xxx-MtKgp2mleNTDL2te3Y-JNBSX2klxVY0
## Source : https://maps.googleapis.com/maps/api/geocode/json?address=Warringah+Mall&key=xxx-MtKgp2mleNTDL2te3Y-JNBSX2klxVY0
## Source : https://maps.googleapis.com/maps/api/geocode/json?address=Westield+Miranda&key=xxx-MtKgp2mleNTDL2te3Y-JNBSX2klxVY0
## Source : https://maps.googleapis.com/maps/api/geocode/json?address=Dee+Why+Grand&key=xxx-MtKgp2mleNTDL2te3Y-JNBSX2klxVY0
## "Dee Why Grand" not uniquely geocoded, using "19 pacific parade, dee why nsw 2099, australia"
## Source : https://maps.googleapis.com/maps/api/geocode/json?address=Westfield+Chatswood&key=xxx-MtKgp2mleNTDL2te3Y-JNBSX2klxVY0
## Source : https://maps.googleapis.com/maps/api/geocode/json?address=Stockland+Piccadilly+Shopping+Centre&key=xxx-MtKgp2mleNTDL2te3Y-JNBSX2klxVY0
## Source : https://maps.googleapis.com/maps/api/geocode/json?address=Westfield+Mount+Druitt&key=xxx-MtKgp2mleNTDL2te3Y-JNBSX2klxVY0
## Source : https://maps.googleapis.com/maps/api/geocode/json?address=South+Village,+580+Princes+Highway+Kirrawee+NSW+2232&key=xxx-MtKgp2mleNTDL2te3Y-JNBSX2klxVY0
## Source : https://maps.googleapis.com/maps/api/geocode/json?address=1+Bay+Street,+Ultimo+NSW+2007&key=xxx-MtKgp2mleNTDL2te3Y-JNBSX2klxVY0
## Source : https://maps.googleapis.com/maps/api/geocode/json?address=Strathfield+Plaza&key=xxx-MtKgp2mleNTDL2te3Y-JNBSX2klxVY0
## Source : https://maps.googleapis.com/maps/api/geocode/json?address=Kings+Cross+Centre&key=xxx-MtKgp2mleNTDL2te3Y-JNBSX2klxVY0
## Source : https://maps.googleapis.com/maps/api/geocode/json?address=Auburn+Cebtral&key=xxx-MtKgp2mleNTDL2te3Y-JNBSX2klxVY0
library(tidyverse)
## ── Attaching packages ─────────────────────────────────────── tidyverse 1.3.0 ──
## ✓ tibble 3.0.4 ✓ purrr 0.3.4
## ✓ tidyr 1.1.2 ✓ stringr 1.4.0
## ✓ readr 1.4.0 ✓ forcats 0.5.0
## ── Conflicts ────────────────────────────────────────── tidyverse_conflicts() ──
## x dplyr::filter() masks stats::filter()
## x dplyr::lag() masks stats::lag()
sc2 <- as_tibble(sc) #make geocoded data into its own dataset
print(sc2, n = 35) #check to see that it is correct. sometimes it will put down the address of a similar named place in another location. Check to see if any values are very different, and then change your dataset to be more specific (like a full street address)
## # A tibble: 36 x 9
## SC Suburb Person Population Income Pop IncomeValue lon lat
## <chr> <chr> <chr> <chr> <chr> <dbl> <dbl> <dbl> <dbl>
## 1 Westfield Eas… Eastga… Tarryn LOW HIGH 860 1833 151. -33.9
## 2 Fashion Spree… Liverp… Tarryn HIGH LOW 27804 1089 151. -33.9
## 3 Neeta City Fairfi… Tarryn LOW LOW 18081 975 151. -33.9
## 4 Marrickville … Marric… Tarryn HIGH HIGH 26592 1814 151. -33.9
## 5 Westfield Pen… Penrith Gemma LOW LOW 10468 2559 151. -33.8
## 6 Westpoint Sho… Blackt… Gemma HIGH LOW 47176 1459 151. -33.8
## 7 Randwick Shop… Randwi… Gemma LOW LOW 13295 1142 151. -33.9
## 8 Kellyville Vi… Kellyv… Gemma HIGH HIGH 29986 2169 151. -33.7
## 9 Rouse Hill To… Rouse … Charl… LOW HIGH 7965 2401 151. -33.7
## 10 Hurstville We… Hurstv… Charl… HIGH LOW 29822 1382 151. -34.0
## 11 Revesby Villa… Revesby Charl… LOW LOW 14176 1516 151. -34.0
## 12 1-3 Brady Str… Mosman Charl… HIGH HIGH 28475 2533 151. -33.8
## 13 Birkenhead Po… Drummo… Corne… LOW HIGH 11950 2353 151. -33.9
## 14 Westfield Hor… Hornsby Corne… HIGH LOW 12551 1383 151. -33.7
## 15 Minto Marketp… Minto Corne… LOW LOW 22168 1696 151. -34.0
## 16 Castle Towers… Castle… Corne… HIGH HIGH 39594 2219 151. -33.7
## 17 Westfield Par… Parram… Yejin LOW HIGH 25798 1739 151. -33.8
## 18 Westfield Bon… Bondi … Yejin HIGH LOW 10986 2026 151. -33.9
## 19 Top Ryde City… Top Ry… Yejin LOW LOW 7929 1604 151. -33.8
## 20 Rhodes Waters… Rhodes Yejin HIGH HIGH 13232 1617 151. -33.8
## 21 Westfield Bur… Burwood Domin… LOW LOW 16030 1398 151. -33.9
## 22 DFO Homebush Homebu… Domin… LOW HIGH 7007 1810 151. -33.9
## 23 Bankstown Cen… Bankst… Domin… HIGH LOW 32113 1120 151. -33.9
## 24 Ashfield Mall Ashfie… Domin… HIGH HIGH 23841 1632 151. -33.9
## 25 Macquarie Cen… Macqua… Ella LOW LOW 8144 1631 151. -33.8
## 26 Warringah Mall Brookv… Ella LOW HIGH 3161 1822 151. -33.8
## 27 Westield Mira… Miranda Ella HIGH LOW 29653 1267 151. -34.0
## 28 Dee Why Grand Dee Why Ella HIGH HIGH 21518 1716 151. -33.8
## 29 Westfield Cha… Chatsw… Kiki HIGH HIGH 24913 1881 151. -33.8
## 30 Stockland Pic… Sydney Kiki LOW HIGH 17252 1949 151. -33.9
## 31 Westfield Mou… Mount … Kiki LOW LOW 16726 1268 151. -33.8
## 32 South Village… Suther… Kiki HIGH LOW 20868 1692 151. -34.0
## 33 1 Bay Street,… Ultimo Vivvy LOW LOW 8845 1230 151. -33.9
## 34 Strathfield P… Strath… Vivvy HIGH HIGH 25813 1892 151. -33.9
## 35 Kings Cross C… Potts … Vivvy LOW HIGH 17252 1949 151. -33.9
## # … with 1 more row
str(sc2)
## tibble [36 × 9] (S3: tbl_df/tbl/data.frame)
## $ SC : chr [1:36] "Westfield Eastgardens" "Fashion Spree Orange Grove" "Neeta City" "Marrickville Metro" ...
## $ Suburb : chr [1:36] "Eastgardens" "Liverpool" "Fairfield" "Marrickville" ...
## $ Person : chr [1:36] "Tarryn" "Tarryn" "Tarryn" "Tarryn" ...
## $ Population : chr [1:36] "LOW" "HIGH" "LOW" "HIGH" ...
## $ Income : chr [1:36] "HIGH" "LOW" "LOW" "HIGH" ...
## $ Pop : num [1:36] 860 27804 18081 26592 10468 ...
## $ IncomeValue: num [1:36] 1833 1089 975 1814 2559 ...
## $ lon : num [1:36] 151 151 151 151 151 ...
## $ lat : num [1:36] -33.9 -33.9 -33.9 -33.9 -33.8 ...
mapview
library(ggplot2)
```r
geocode("Parramatta") #check coordinates of where you want your centre of the map to be. I changed from Sydney to Parramatta to get all of the values in.
## Source : https://maps.googleapis.com/maps/api/geocode/json?address=Parramatta&key=xxx-MtKgp2mleNTDL2te3Y-JNBSX2klxVY0
## # A tibble: 1 x 2
## lon lat
## <dbl> <dbl>
## 1 151. -33.8
map <- get_googlemap(center = c(151.0,-33.8), zoom = 10) #using the results from the code above, sub in your coordinates to set the centre of the map.
## Source : https://maps.googleapis.com/maps/api/staticmap?center=-33.8,151&zoom=10&size=640x640&scale=2&maptype=terrain&key=xxx-MtKgp2mleNTDL2te3Y-JNBSX2klxVY0
# ggmap has the same syntax as ggplot
ggmap(map) +
geom_point(data = sc2, aes(x= lon, y = lat, size = Pop, color = IncomeValue)) #this line adds plot points. Set a variable to "size" and "colour" to map the values of that variable using size and colour of plots. Make sure this is put inside the aes() brackets.
ggmap(map) +
geom_point(data = sc2, aes(x= lon, y = lat, size = Pop, color = IncomeValue))
ggmap(map) +
geom_point(data = sc2, aes(x= lon, y = lat, shape = Population, color = Income), size = 3) +
labs(caption = "(based on data collected by the Australian Bureau of Statistics in the 2016 census)", x = "Longitude", y = "Latitude")
str(sc2)
## tibble [36 × 9] (S3: tbl_df/tbl/data.frame)
## $ SC : chr [1:36] "Westfield Eastgardens" "Fashion Spree Orange Grove" "Neeta City" "Marrickville Metro" ...
## $ Suburb : chr [1:36] "Eastgardens" "Liverpool" "Fairfield" "Marrickville" ...
## $ Person : chr [1:36] "Tarryn" "Tarryn" "Tarryn" "Tarryn" ...
## $ Population : chr [1:36] "LOW" "HIGH" "LOW" "HIGH" ...
## $ Income : chr [1:36] "HIGH" "LOW" "LOW" "HIGH" ...
## $ Pop : num [1:36] 860 27804 18081 26592 10468 ...
## $ IncomeValue: num [1:36] 1833 1089 975 1814 2559 ...
## $ lon : num [1:36] 151 151 151 151 151 ...
## $ lat : num [1:36] -33.9 -33.9 -33.9 -33.9 -33.8 ...
mapview
library(ggplot2)
geocode("Parramatta") #check coordinates of where you want your centre of the map to be. I changed from Sydney to Parramatta to get all of the values in.
map <- get_googlemap(center = c(151.0,-33.8), zoom = 10) #using the results from the code above, sub in your coordinates to set the centre of the map.
# ggmap has the same syntax as ggplot
ggmap(map) +
geom_point(data = sc2, aes(x= lon, y = lat, size = Pop, color = IncomeValue)) #this line adds plot points. Set a variable to "size" and "colour" to map the values of that variable using size and colour of plots. Make sure this is put inside the aes() brackets.
ggmap(map) +
geom_point(data = sc2, aes(x= lon, y = lat, size = Pop, color = IncomeValue))
ggmap(map) +
geom_point(data = sc2, aes(x= lon, y = lat, shape = Population, color = Income), size = 3) +
labs(caption = "(based on data collected by the Australian Bureau of Statistics in the 2016 census)", x = "Longitude", y = "Latitude")
library(readr)
hab1 <- read_csv("withfactors.csv", col_types = cols(Population = col_factor(levels = c("HIGH",
"LOW")), Income = col_factor(levels = c("HIGH",
"LOW"))))
## Warning: Missing column names filled in: 'X16' [16], 'X17' [17], 'X18' [18],
## 'X19' [19]
library(REdaS)
## Loading required package: grid
cormatrix <- cor(hab1[1:36,4:15])
cormatrix
## carparks escalators elevators levels toilets
## carparks 1.00000000 0.1104018 -0.03651974 0.27488040 0.36134990
## escalators 0.11040178 1.0000000 0.73954196 0.50882483 0.71726211
## elevators -0.03651974 0.7395420 1.00000000 0.44497674 0.47146012
## levels 0.27488040 0.5088248 0.44497674 1.00000000 0.38411470
## toilets 0.36134990 0.7172621 0.47146012 0.38411470 1.00000000
## hours 0.12710142 0.1515900 0.09310969 0.01157244 0.08566331
## eateries 0.40250032 0.4862020 0.31493815 0.49847231 0.67115063
## trainstations 0.18535517 0.3602172 0.22379553 0.27026228 0.34410700
## busstops 0.12063463 0.1247951 0.15460860 -0.08987194 -0.02223400
## cinema 0.44525146 0.3371226 0.22596641 0.68586483 0.50876630
## greenspace -0.05889480 -0.2018897 -0.07751053 -0.22812489 -0.14558848
## foodcourts 0.42506151 0.4548239 0.07409223 0.31467775 0.37917306
## hours eateries trainstations busstops cinema
## carparks 0.12710142 0.4025003 0.1853552 0.12063463 0.44525146
## escalators 0.15159000 0.4862020 0.3602172 0.12479512 0.33712259
## elevators 0.09310969 0.3149381 0.2237955 0.15460860 0.22596641
## levels 0.01157244 0.4984723 0.2702623 -0.08987194 0.68586483
## toilets 0.08566331 0.6711506 0.3441070 -0.02223400 0.50876630
## hours 1.00000000 0.2923818 0.1756103 0.39344563 0.09786211
## eateries 0.29238185 1.0000000 0.2516909 0.11067248 0.69984069
## trainstations 0.17561033 0.2516909 1.0000000 -0.12707898 0.30303450
## busstops 0.39344563 0.1106725 -0.1270790 1.00000000 -0.14285714
## cinema 0.09786211 0.6998407 0.3030345 -0.14285714 1.00000000
## greenspace -0.28375623 -0.3502501 -0.2115139 -0.17131295 -0.14151940
## foodcourts 0.10887017 0.4216912 0.3893209 0.06992765 0.34963827
## greenspace foodcourts
## carparks -0.05889480 0.42506151
## escalators -0.20188969 0.45482390
## elevators -0.07751053 0.07409223
## levels -0.22812489 0.31467775
## toilets -0.14558848 0.37917306
## hours -0.28375623 0.10887017
## eateries -0.35025006 0.42169121
## trainstations -0.21151389 0.38932091
## busstops -0.17131295 0.06992765
## cinema -0.14151940 0.34963827
## greenspace 1.00000000 0.06131806
## foodcourts 0.06131806 1.00000000
library(ggcorrplot)
ggcorrplot(cormatrix)
ggcorrplot_clustered <- ggcorrplot(cormatrix, hc.order = TRUE, type ="lower")
ggcorrplot_clustered
Hierarchical clustering of correlation. Eelvators and escalators highly correlated, eateries and toilets, cinema and eateries and toilets and escaltors are all highly correlated.
bart_spher(cormatrix)
## Bartlett's Test of Sphericity
##
## Call: bart_spher(x = cormatrix)
##
## X2 = 325.279
## df = 66
## p-value < 2.22e-16
Bartlett’s Test of Sphericity = testing assumptions to see if the data is correlated.
p-value = <2.222e<16
Make PCs
pca1 <- prcomp(cor(hab1[1:36,4:15]), scale = TRUE)
pca2 <- princomp(cor(hab1[1:36,4:15]), cor = TRUE)
summary(pca1)
## Importance of components:
## PC1 PC2 PC3 PC4 PC5 PC6 PC7
## Standard deviation 2.3153 1.4838 1.3814 0.96836 0.83908 0.59568 0.55923
## Proportion of Variance 0.4467 0.1835 0.1590 0.07814 0.05867 0.02957 0.02606
## Cumulative Proportion 0.4467 0.6302 0.7892 0.86734 0.92601 0.95558 0.98164
## PC8 PC9 PC10 PC11 PC12
## Standard deviation 0.3811 0.21125 0.15612 0.07828 6.375e-17
## Proportion of Variance 0.0121 0.00372 0.00203 0.00051 0.000e+00
## Cumulative Proportion 0.9937 0.99746 0.99949 1.00000 1.000e+00
PC1 explains 44.67% of the variance, PC2 explains 18.35% of the total variance, PC3 explains 15.90% and PC4 explains 7.81% of the total variance.These 4 PCs explain a cumulative proportion of 86.7% of the variance.
screeplot(pca1, type = "lines")
pca1$x
## PC1 PC2 PC3 PC4 PC5
## carparks 0.1359197 1.42762788 -2.12064115 -0.6480406 0.42151020
## escalators -1.7758270 -0.98825800 1.47672596 0.4624595 0.84534876
## elevators -0.3488096 -1.45542020 2.42897906 -0.3110332 0.01612957
## levels -1.6302632 0.45123704 0.63126664 -0.8687229 -0.82116336
## toilets -1.9955015 0.02737616 0.53731394 -0.1275671 0.62717405
## hours 1.9504040 -2.18421583 -1.54202171 0.4378214 -0.98924196
## eateries -2.0217986 -0.31036094 -0.88292413 -0.7082315 0.02691283
## trainstations -0.4321372 0.40809526 0.08074933 2.2304057 -1.21981355
## busstops 3.3206029 -2.24103373 -0.84759577 -0.4980453 0.75193287
## cinema -1.7808635 1.08700961 -0.49980569 -1.0038613 -0.92467426
## greenspace 5.1417437 2.39232169 1.65758265 -0.3054750 -0.10713926
## foodcourts -0.5634696 1.38562107 -0.91962911 1.3402903 1.37302411
## PC6 PC7 PC8 PC9 PC10
## carparks -0.24735751 1.07550188 -0.52082425 -0.06842893 0.0307474622
## escalators 0.02636193 -0.09866588 -0.25190649 0.16329447 -0.1002073999
## elevators 0.05055560 0.38441168 -0.29792874 -0.40997153 0.0003711169
## levels 1.13341732 0.01701822 -0.17031582 0.30597723 0.1837423642
## toilets -1.15276440 0.07944638 0.16885678 0.31872145 -0.0163579956
## hours -0.37532273 -0.70890312 -0.56631294 0.06364605 -0.0118396531
## eateries -0.40709415 -0.55098661 0.46695602 -0.24158190 0.2843990468
## trainstations -0.04749158 0.67870695 0.37529164 -0.01684887 0.0341781632
## busstops 0.55132449 0.44108681 0.60034682 0.08590062 -0.0629123851
## cinema 0.10271361 -0.29838173 0.25561152 -0.10693748 -0.3663499310
## greenspace -0.33440181 -0.35220951 0.03814829 0.02502054 0.0469579081
## foodcourts 0.70005923 -0.66702509 -0.09792283 -0.11879164 -0.0227286967
## PC11 PC12
## carparks 0.036258646 3.191891e-16
## escalators 0.178261906 -5.689893e-16
## elevators -0.080050390 -1.665335e-16
## levels -0.036670109 5.342948e-16
## toilets -0.120354080 1.838807e-16
## hours -0.018750742 3.330669e-16
## eateries 0.068548932 -6.383782e-16
## trainstations 0.017061656 -1.040834e-17
## busstops -0.008841529 -5.273559e-16
## cinema 0.006878692 3.903128e-17
## greenspace 0.032259949 1.387779e-16
## foodcourts -0.074602930 3.330669e-16
habpcscores <- cbind(hab1[1:36,1:15], pca1$x)
## Warning in data.frame(..., check.names = FALSE): row names were found from a
## short variable and have been discarded
habpcscores
## Shopping Centre Population Income carparks escalators
## 1 Westfield Eastgardens LOW HIGH 1 9
## 2 Fashion Spree Orange Grove HIGH LOW 3 0
## 3 Neeta City LOW LOW 1 10
## 4 Marrickville Metro HIGH HIGH 1 8
## 5 Kellyville Village LOW HIGH 1 1
## 6 Blacktown Westpoint HIGH LOW 5 9
## 7 Penrith Westfield LOW LOW 6 19
## 8 Royal Randwick Shopping Centre HIGH HIGH 1 0
## 9 Rouse Hill Town Centre LOW HIGH 4 5
## 10 Hurstville Westfield HIGH LOW 10 14
## 11 Revesby Village Centre LOW LOW 1 6
## 12 Bridgepoint Shopping Centre HIGH HIGH 1 2
## 13 Birkenhead Point LOW HIGH 3 16
## 14 Westfield Shopping Centre HIGH LOW 1 21
## 15 Minto Marketplace LOW LOW 2 2
## 16 Castle Towers HIGH HIGH 3 1
## 17 Westfield Bondi Junction LOW HIGH 3 49
## 18 Rhodes Waterside HIGH LOW 1 3
## 19 Top Ryde LOW LOW 1 6
## 20 Westfield Parramatta HIGH HIGH 2 4
## 21 DFO Homebush LOW HIGH 2 8
## 22 Bankstown Central HIGH LOW 2 24
## 23 Westfield Burwood LOW LOW 5 21
## 24 Ashfield Mall HIGH HIGH 1 6
## 25 Macquarie Centre LOW LOW 5 17
## 26 Warringah Mall LOW HIGH 6 7
## 27 Westield Miranda HIGH LOW 3 28
## 28 Dee Why Grand HIGH HIGH 1 4
## 29 Stockland Piccadilly Shopping Centre LOW HIGH 1 2
## 30 South Village HIGH LOW 1 12
## 31 Westfield Mount Druitt LOW LOW 7 4
## 32 Westfield Chatswood HIGH HIGH 1 64
## 33 Kings Cross Centre LOW HIGH 1 0
## 34 Auburn Central HIGH LOW 1 6
## 35 Broadway Shopping Centre LOW LOW 2 2
## 36 Strathfield Plaza HIGH HIGH 1 1
## elevators levels toilets hours eateries trainstations busstops cinema
## 1 6 3 21 12.0 52 1 1 1
## 2 0 1 6 11.0 6 0 1 0
## 3 3 3 5 12.0 19 1 1 0
## 4 2 2 2 12.0 20 0 1 0
## 5 2 1 1 12.0 12 0 1 0
## 6 7 9 6 12.0 32 1 1 1
## 7 4 5 8 12.0 45 1 1 1
## 8 2 2 6 12.0 19 0 1 0
## 9 2 4 10 12.0 31 1 1 1
## 10 3 6 14 11.5 48 1 1 1
## 11 1 3 4 12.0 3 1 1 0
## 12 2 2 2 12.0 6 0 1 0
## 13 5 4 6 9.5 18 0 1 0
## 14 4 4 5 12.0 65 1 1 1
## 15 3 1 4 8.5 16 1 1 0
## 16 1 4 7 12.0 61 0 1 1
## 17 4 7 24 11.5 47 1 1 1
## 18 1 7 3 10.0 38 0 1 1
## 19 6 7 6 12.0 63 0 1 1
## 20 3 6 3 9.5 26 1 1 1
## 21 3 2 12 8.0 15 0 1 0
## 22 5 3 17 12.0 46 1 1 1
## 23 5 5 23 11.5 41 1 1 1
## 24 6 4 7 12.0 18 1 1 0
## 25 5 5 24 11.5 104 1 1 1
## 26 6 3 13 12.0 44 0 1 1
## 27 8 5 22 11.5 74 1 1 1
## 28 2 2 6 11.0 8 0 1 0
## 29 4 2 3 11.0 10 1 1 0
## 30 10 3 4 11.5 11 1 1 0
## 31 3 3 9 12.0 47 1 1 1
## 32 29 8 24 11.5 57 1 1 1
## 33 2 3 1 12.0 16 1 1 1
## 34 1 3 3 12.0 14 1 1 0
## 35 6 5 5 9.0 7 0 1 1
## 36 0 5 10 8.5 17 1 0 1
## greenspace foodcourts PC1 PC2 PC3 PC4
## 1 0 1 0.1359197 1.42762788 -2.12064115 -0.6480406
## 2 0 1 -1.7758270 -0.98825800 1.47672596 0.4624595
## 3 0 1 -0.3488096 -1.45542020 2.42897906 -0.3110332
## 4 0 1 -1.6302632 0.45123704 0.63126664 -0.8687229
## 5 0 0 -1.9955015 0.02737616 0.53731394 -0.1275671
## 6 0 2 1.9504040 -2.18421583 -1.54202171 0.4378214
## 7 0 2 -2.0217986 -0.31036094 -0.88292413 -0.7082315
## 8 0 0 -0.4321372 0.40809526 0.08074933 2.2304057
## 9 1 2 3.3206029 -2.24103373 -0.84759577 -0.4980453
## 10 0 2 -1.7808635 1.08700961 -0.49980569 -1.0038613
## 11 0 1 5.1417437 2.39232169 1.65758265 -0.3054750
## 12 3 1 -0.5634696 1.38562107 -0.91962911 1.3402903
## 13 1 2 0.1359197 1.42762788 -2.12064115 -0.6480406
## 14 0 3 -1.7758270 -0.98825800 1.47672596 0.4624595
## 15 1 2 -0.3488096 -1.45542020 2.42897906 -0.3110332
## 16 0 1 -1.6302632 0.45123704 0.63126664 -0.8687229
## 17 0 3 -1.9955015 0.02737616 0.53731394 -0.1275671
## 18 0 1 1.9504040 -2.18421583 -1.54202171 0.4378214
## 19 0 1 -2.0217986 -0.31036094 -0.88292413 -0.7082315
## 20 0 1 -0.4321372 0.40809526 0.08074933 2.2304057
## 21 1 1 3.3206029 -2.24103373 -0.84759577 -0.4980453
## 22 0 2 -1.7808635 1.08700961 -0.49980569 -1.0038613
## 23 1 2 5.1417437 2.39232169 1.65758265 -0.3054750
## 24 0 1 -0.5634696 1.38562107 -0.91962911 1.3402903
## 25 0 2 0.1359197 1.42762788 -2.12064115 -0.6480406
## 26 1 2 -1.7758270 -0.98825800 1.47672596 0.4624595
## 27 0 1 -0.3488096 -1.45542020 2.42897906 -0.3110332
## 28 0 1 -1.6302632 0.45123704 0.63126664 -0.8687229
## 29 0 0 -1.9955015 0.02737616 0.53731394 -0.1275671
## 30 1 2 1.9504040 -2.18421583 -1.54202171 0.4378214
## 31 0 1 -2.0217986 -0.31036094 -0.88292413 -0.7082315
## 32 0 1 -0.4321372 0.40809526 0.08074933 2.2304057
## 33 1 1 3.3206029 -2.24103373 -0.84759577 -0.4980453
## 34 0 1 -1.7808635 1.08700961 -0.49980569 -1.0038613
## 35 1 0 5.1417437 2.39232169 1.65758265 -0.3054750
## 36 1 1 -0.5634696 1.38562107 -0.91962911 1.3402903
## PC5 PC6 PC7 PC8 PC9 PC10
## 1 0.42151020 -0.24735751 1.07550188 -0.52082425 -0.06842893 0.0307474622
## 2 0.84534876 0.02636193 -0.09866588 -0.25190649 0.16329447 -0.1002073999
## 3 0.01612957 0.05055560 0.38441168 -0.29792874 -0.40997153 0.0003711169
## 4 -0.82116336 1.13341732 0.01701822 -0.17031582 0.30597723 0.1837423642
## 5 0.62717405 -1.15276440 0.07944638 0.16885678 0.31872145 -0.0163579956
## 6 -0.98924196 -0.37532273 -0.70890312 -0.56631294 0.06364605 -0.0118396531
## 7 0.02691283 -0.40709415 -0.55098661 0.46695602 -0.24158190 0.2843990468
## 8 -1.21981355 -0.04749158 0.67870695 0.37529164 -0.01684887 0.0341781632
## 9 0.75193287 0.55132449 0.44108681 0.60034682 0.08590062 -0.0629123851
## 10 -0.92467426 0.10271361 -0.29838173 0.25561152 -0.10693748 -0.3663499310
## 11 -0.10713926 -0.33440181 -0.35220951 0.03814829 0.02502054 0.0469579081
## 12 1.37302411 0.70005923 -0.66702509 -0.09792283 -0.11879164 -0.0227286967
## 13 0.42151020 -0.24735751 1.07550188 -0.52082425 -0.06842893 0.0307474622
## 14 0.84534876 0.02636193 -0.09866588 -0.25190649 0.16329447 -0.1002073999
## 15 0.01612957 0.05055560 0.38441168 -0.29792874 -0.40997153 0.0003711169
## 16 -0.82116336 1.13341732 0.01701822 -0.17031582 0.30597723 0.1837423642
## 17 0.62717405 -1.15276440 0.07944638 0.16885678 0.31872145 -0.0163579956
## 18 -0.98924196 -0.37532273 -0.70890312 -0.56631294 0.06364605 -0.0118396531
## 19 0.02691283 -0.40709415 -0.55098661 0.46695602 -0.24158190 0.2843990468
## 20 -1.21981355 -0.04749158 0.67870695 0.37529164 -0.01684887 0.0341781632
## 21 0.75193287 0.55132449 0.44108681 0.60034682 0.08590062 -0.0629123851
## 22 -0.92467426 0.10271361 -0.29838173 0.25561152 -0.10693748 -0.3663499310
## 23 -0.10713926 -0.33440181 -0.35220951 0.03814829 0.02502054 0.0469579081
## 24 1.37302411 0.70005923 -0.66702509 -0.09792283 -0.11879164 -0.0227286967
## 25 0.42151020 -0.24735751 1.07550188 -0.52082425 -0.06842893 0.0307474622
## 26 0.84534876 0.02636193 -0.09866588 -0.25190649 0.16329447 -0.1002073999
## 27 0.01612957 0.05055560 0.38441168 -0.29792874 -0.40997153 0.0003711169
## 28 -0.82116336 1.13341732 0.01701822 -0.17031582 0.30597723 0.1837423642
## 29 0.62717405 -1.15276440 0.07944638 0.16885678 0.31872145 -0.0163579956
## 30 -0.98924196 -0.37532273 -0.70890312 -0.56631294 0.06364605 -0.0118396531
## 31 0.02691283 -0.40709415 -0.55098661 0.46695602 -0.24158190 0.2843990468
## 32 -1.21981355 -0.04749158 0.67870695 0.37529164 -0.01684887 0.0341781632
## 33 0.75193287 0.55132449 0.44108681 0.60034682 0.08590062 -0.0629123851
## 34 -0.92467426 0.10271361 -0.29838173 0.25561152 -0.10693748 -0.3663499310
## 35 -0.10713926 -0.33440181 -0.35220951 0.03814829 0.02502054 0.0469579081
## 36 1.37302411 0.70005923 -0.66702509 -0.09792283 -0.11879164 -0.0227286967
## PC11 PC12
## 1 0.036258646 3.191891e-16
## 2 0.178261906 -5.689893e-16
## 3 -0.080050390 -1.665335e-16
## 4 -0.036670109 5.342948e-16
## 5 -0.120354080 1.838807e-16
## 6 -0.018750742 3.330669e-16
## 7 0.068548932 -6.383782e-16
## 8 0.017061656 -1.040834e-17
## 9 -0.008841529 -5.273559e-16
## 10 0.006878692 3.903128e-17
## 11 0.032259949 1.387779e-16
## 12 -0.074602930 3.330669e-16
## 13 0.036258646 3.191891e-16
## 14 0.178261906 -5.689893e-16
## 15 -0.080050390 -1.665335e-16
## 16 -0.036670109 5.342948e-16
## 17 -0.120354080 1.838807e-16
## 18 -0.018750742 3.330669e-16
## 19 0.068548932 -6.383782e-16
## 20 0.017061656 -1.040834e-17
## 21 -0.008841529 -5.273559e-16
## 22 0.006878692 3.903128e-17
## 23 0.032259949 1.387779e-16
## 24 -0.074602930 3.330669e-16
## 25 0.036258646 3.191891e-16
## 26 0.178261906 -5.689893e-16
## 27 -0.080050390 -1.665335e-16
## 28 -0.036670109 5.342948e-16
## 29 -0.120354080 1.838807e-16
## 30 -0.018750742 3.330669e-16
## 31 0.068548932 -6.383782e-16
## 32 0.017061656 -1.040834e-17
## 33 -0.008841529 -5.273559e-16
## 34 0.006878692 3.903128e-17
## 35 0.032259949 1.387779e-16
## 36 -0.074602930 3.330669e-16
pc1aov <- aov(habpcscores$PC1 ~ busstops*greenspace*hours, data = habpcscores)
summary(pc1aov)
## Df Sum Sq Mean Sq F value Pr(>F)
## busstops 1 0.33 0.327 0.072 0.7906
## greenspace 1 24.18 24.179 5.309 0.0281 *
## hours 1 8.71 8.707 1.912 0.1766
## greenspace:hours 1 2.51 2.509 0.551 0.4635
## Residuals 31 141.17 4.554
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
pc2aov <- aov(habpcscores$PC2 ~ carparks*greenspace*cinema, data = habpcscores)
summary(pc2aov)
## Df Sum Sq Mean Sq F value Pr(>F)
## carparks 1 0.00 0.0008 0.000 0.986
## greenspace 1 0.01 0.0118 0.005 0.946
## cinema 1 0.08 0.0750 0.030 0.864
## carparks:greenspace 1 0.02 0.0210 0.008 0.928
## carparks:cinema 1 1.46 1.4556 0.579 0.453
## greenspace:cinema 1 0.25 0.2531 0.101 0.753
## carparks:greenspace:cinema 1 0.44 0.4435 0.176 0.678
## Residuals 28 70.39 2.5140
pc3aov <- aov(habpcscores$PC3 ~ escalators*elevators, data = habpcscores)
summary(pc3aov)
## Df Sum Sq Mean Sq F value Pr(>F)
## escalators 1 0.19 0.1935 0.102 0.751
## elevators 1 1.12 1.1213 0.591 0.448
## escalators:elevators 1 0.98 0.9845 0.519 0.476
## Residuals 32 60.67 1.8960
pc4aov <- aov(habpcscores$PC4 ~ foodcourts*trainstations*hours, data = habpcscores)
summary(pc4aov)
## Df Sum Sq Mean Sq F value Pr(>F)
## foodcourts 1 1.168 1.168 1.414 0.2444
## trainstations 1 0.195 0.195 0.236 0.6307
## hours 1 0.469 0.469 0.568 0.4575
## foodcourts:trainstations 1 0.127 0.127 0.154 0.6978
## foodcourts:hours 1 0.003 0.003 0.003 0.9545
## trainstations:hours 1 3.892 3.892 4.714 0.0385 *
## foodcourts:trainstations:hours 1 1.973 1.973 2.390 0.1333
## Residuals 28 23.118 0.826
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
In PC1, greenspace is significant (p = 0.0281).
In PC4, we see an interaction between presence of train stations and opening hours on Thursdays with a p-value of 0.0385.
library(FactoMineR)
sc_pca <- PCA(hab1 [1:36,4:15])
summary(sc_pca)
##
## Call:
## PCA(X = hab1[1:36, 4:15])
##
##
## Eigenvalues
## Dim.1 Dim.2 Dim.3 Dim.4 Dim.5 Dim.6 Dim.7
## Variance 4.351 1.583 1.422 1.087 0.958 0.650 0.618
## % of var. 36.262 13.190 11.852 9.057 7.986 5.416 5.151
## Cumulative % of var. 36.262 49.452 61.304 70.360 78.347 83.763 88.913
## Dim.8 Dim.9 Dim.10 Dim.11 Dim.12
## Variance 0.520 0.366 0.214 0.150 0.080
## % of var. 4.337 3.049 1.784 1.249 0.668
## Cumulative % of var. 93.251 96.299 98.084 99.332 100.000
##
## Individuals (the 10 first)
## Dist Dim.1 ctr cos2 Dim.2 ctr cos2 Dim.3
## 1 | 2.573 | 1.242 0.985 0.233 | 0.601 0.633 0.054 | -0.274
## 2 | 2.991 | -2.424 3.751 0.657 | 0.661 0.767 0.049 | 0.581
## 3 | 2.092 | -1.045 0.698 0.250 | 1.084 2.062 0.269 | -0.335
## 4 | 2.668 | -1.992 2.533 0.558 | 1.303 2.981 0.239 | -0.096
## 5 | 3.548 | -2.932 5.487 0.683 | 1.506 3.980 0.180 | -0.242
## 6 | 3.358 | 2.011 2.582 0.359 | -0.148 0.038 0.002 | 0.796
## 7 | 2.562 | 1.863 2.217 0.529 | 0.184 0.059 0.005 | 1.356
## 8 | 3.154 | -2.400 3.678 0.579 | 1.378 3.332 0.191 | -0.287
## 9 | 2.138 | 0.589 0.222 0.076 | -0.443 0.344 0.043 | 1.343
## 10 | 4.112 | 2.589 4.277 0.396 | -0.405 0.289 0.010 | 2.419
## ctr cos2
## 1 0.147 0.011 |
## 2 0.659 0.038 |
## 3 0.219 0.026 |
## 4 0.018 0.001 |
## 5 0.115 0.005 |
## 6 1.237 0.056 |
## 7 3.593 0.280 |
## 8 0.161 0.008 |
## 9 3.524 0.395 |
## 10 11.429 0.346 |
##
## Variables (the 10 first)
## Dim.1 ctr cos2 Dim.2 ctr cos2 Dim.3 ctr
## carparks | 0.484 5.385 0.234 | -0.101 0.641 0.010 | 0.641 28.847
## escalators | 0.785 14.173 0.617 | 0.103 0.674 0.011 | -0.472 15.633
## elevators | 0.569 7.453 0.324 | 0.152 1.453 0.023 | -0.693 33.766
## levels | 0.718 11.848 0.516 | -0.229 3.325 0.053 | -0.104 0.763
## toilets | 0.800 14.698 0.640 | -0.089 0.505 0.008 | -0.138 1.348
## hours | 0.269 1.662 0.072 | 0.711 31.918 0.505 | 0.259 4.710
## eateries | 0.820 15.468 0.673 | 0.081 0.414 0.007 | 0.199 2.787
## trainstations | 0.517 6.149 0.268 | -0.080 0.406 0.006 | 0.025 0.043
## busstops | 0.084 0.163 0.007 | 0.788 39.209 0.621 | 0.099 0.684
## cinema | 0.747 12.818 0.558 | -0.302 5.771 0.091 | 0.246 4.263
## cos2
## carparks 0.410 |
## escalators 0.222 |
## elevators 0.480 |
## levels 0.011 |
## toilets 0.019 |
## hours 0.067 |
## eateries 0.040 |
## trainstations 0.001 |
## busstops 0.010 |
## cinema 0.061 |
sc_pca$eig
## eigenvalue percentage of variance cumulative percentage of variance
## comp 1 4.35145114 36.2620928 36.26209
## comp 2 1.58280047 13.1900039 49.45210
## comp 3 1.42218382 11.8515319 61.30363
## comp 4 1.08682219 9.0568516 70.36048
## comp 5 0.95832816 7.9860680 78.34655
## comp 6 0.64994119 5.4161765 83.76272
## comp 7 0.61808282 5.1506902 88.91341
## comp 8 0.52045800 4.3371500 93.25056
## comp 9 0.36583013 3.0485844 96.29915
## comp 10 0.21413594 1.7844662 98.08362
## comp 11 0.14983211 1.2486009 99.33222
## comp 12 0.08013402 0.6677835 100.00000
fa1 <- factanal(hab1[1:36,4:15], 4, rotation = "varimax")
fa1
##
## Call:
## factanal(x = hab1[1:36, 4:15], factors = 4, rotation = "varimax")
##
## Uniquenesses:
## carparks escalators elevators levels toilets
## 0.651 0.005 0.357 0.394 0.346
## hours eateries trainstations busstops cinema
## 0.791 0.151 0.785 0.795 0.070
## greenspace foodcourts
## 0.760 0.005
##
## Loadings:
## Factor1 Factor2 Factor3 Factor4
## carparks 0.438 0.357 0.151
## escalators 0.924 0.114 0.297 0.199
## elevators 0.788 0.104
## levels 0.454 0.606 0.144 -0.113
## toilets 0.587 0.394 0.208 0.334
## hours 0.444
## eateries 0.254 0.672 0.236 0.527
## trainstations 0.259 0.212 0.321
## busstops -0.172 0.409
## cinema 0.206 0.930 0.150
## greenspace -0.153 -0.163 0.141 -0.412
## foodcourts 0.146 0.189 0.967
##
## Factor1 Factor2 Factor3 Factor4
## SS loadings 2.254 2.193 1.430 1.014
## Proportion Var 0.188 0.183 0.119 0.084
## Cumulative Var 0.188 0.371 0.490 0.574
##
## Test of the hypothesis that 4 factors are sufficient.
## The chi square statistic is 26.01 on 24 degrees of freedom.
## The p-value is 0.352
library(vegan)
## Loading required package: permute
## Loading required package: lattice
## This is vegan 2.5-6
library(MASS)
##
## Attaching package: 'MASS'
## The following object is masked from 'package:dplyr':
##
## select
library(ggplot2)
library(readr)
full <- read_csv("full.csv", col_types = cols(Population = col_factor(levels = c("LOW",
"HIGH")), Income = col_factor(levels = c("LOW",
"HIGH")), PopulationValue = col_number(),
IncomeValue = col_number()))
library(readr)
assemblage <- read_csv("assemblage.csv")
## Warning: Missing column names filled in: 'X38' [38], 'X39' [39]
##
## ── Column specification ────────────────────────────────────────────────────────
## cols(
## .default = col_double(),
## Species = col_character(),
## X38 = col_logical(),
## X39 = col_logical()
## )
## ℹ Use `spec()` for the full column specifications.
library(readr)
species <- read_csv("species.csv")
##
## ── Column specification ────────────────────────────────────────────────────────
## cols(
## .default = col_double(),
## ShoppingCentre = col_character(),
## Suburb = col_character(),
## Population = col_character(),
## Income = col_character(),
## PopulationValue = col_number(),
## IncomeValue = col_number()
## )
## ℹ Use `spec()` for the full column specifications.
spp1 <- (na.omit( species
[2:32,4:20])) #square brackets used to select only numeric data
library(readr)
spc <- read_csv("spp.csv")
## Warning: Missing column names filled in: 'X19' [19], 'X20' [20]
##
## ── Column specification ────────────────────────────────────────────────────────
## cols(
## .default = col_double(),
## `Shopping Centre` = col_character(),
## Suburb = col_character(),
## Population = col_character(),
## Income = col_character(),
## X19 = col_logical(),
## X20 = col_logical()
## )
## ℹ Use `spec()` for the full column specifications.
spp <- spc[1:36,5:18]
library(readr)
factors <- read_csv("factors.csv", col_types = cols(Population = col_factor(levels = c("LOW",
"HIGH")), Income = col_factor(levels = c("LOW",
"HIGH"))))
## Warning: Missing column names filled in: 'X1' [1], 'X7' [7], 'X8' [8], 'X9' [9],
## 'X10' [10], 'X11' [11], 'X12' [12], 'X13' [13], 'X14' [14], 'X15' [15],
## 'X16' [16], 'X17' [17], 'X18' [18], 'X19' [19], 'X20' [20], 'X21' [21],
## 'X22' [22], 'X23' [23], 'X24' [24], 'X25' [25], 'X26' [26], 'X27' [27],
## 'X28' [28], 'X29' [29], 'X30' [30], 'X31' [31], 'X32' [32], 'X33' [33],
## 'X34' [34], 'X35' [35], 'X36' [36], 'X37' [37], 'X38' [38], 'X39' [39]
str(spp)
## tibble [36 × 14] (S3: tbl_df/tbl/data.frame)
## $ Department: num [1:36] 4 0 0 2 0 5 3 0 2 2 ...
## $ Grocery : num [1:36] 5 0 4 6 3 10 3 3 5 7 ...
## $ Accessory : num [1:36] 22 14 6 5 1 23 30 3 22 21 ...
## $ Clothing : num [1:36] 36 31 15 7 1 37 55 9 34 30 ...
## $ Pharmacy : num [1:36] 4 0 1 1 1 2 3 1 2 3 ...
## $ Pet : num [1:36] 1 0 0 0 0 0 2 0 0 0 ...
## $ Toy : num [1:36] 2 0 0 1 0 1 1 0 3 1 ...
## $ Specialty : num [1:36] 7 0 1 2 4 6 25 2 2 1 ...
## $ Thrift : num [1:36] 0 0 0 0 0 0 0 0 0 0 ...
## $ Services : num [1:36] 63 2 26 32 18 63 62 22 36 44 ...
## $ Homeware : num [1:36] 6 5 0 1 0 6 4 1 8 3 ...
## $ Technology: num [1:36] 8 1 4 6 1 10 10 1 9 9 ...
## $ Cosmetic : num [1:36] 2 0 0 0 0 5 4 0 5 8 ...
## $ Fresh Food: num [1:36] 52 6 19 20 15 53 70 19 6 43 ...
TransSpp <- spp^(1/4)
spp.dis <- vegdist(TransSpp, method ="bray")
spp.mds0 <- isoMDS(spp.dis)
## initial value 13.795483
## iter 5 value 10.457427
## iter 10 value 10.286952
## iter 15 value 9.983493
## final value 9.873226
## converged
Plotting nMDS
plot <- ordiplot(spp.mds0, type = "t")
## species scores not available
pchs <- c(0:5)
Pop_Factor <- factors$Population
nMDS_spp<- metaMDS(TransSpp, distance="bray", k=2)
## Run 0 stress 0.0969794
## Run 1 stress 0.126783
## Run 2 stress 0.1072847
## Run 3 stress 0.1071067
## Run 4 stress 0.105539
## Run 5 stress 0.1250816
## Run 6 stress 0.1342471
## Run 7 stress 0.1001295
## Run 8 stress 0.1066917
## Run 9 stress 0.1062038
## Run 10 stress 0.1062039
## Run 11 stress 0.1066749
## Run 12 stress 0.1357998
## Run 13 stress 0.1199346
## Run 14 stress 0.1198358
## Run 15 stress 0.1274054
## Run 16 stress 0.1171052
## Run 17 stress 0.1166742
## Run 18 stress 0.1072848
## Run 19 stress 0.1062039
## Run 20 stress 0.1267099
## *** No convergence -- monoMDS stopping criteria:
## 20: stress ratio > sratmax
plot_spp <- ordiplot(nMDS_spp, display = "sites", type = "n")
points(nMDS_spp, col = "black", pch = pchs[Pop_Factor])
legend("topright", bty = "n", legend = levels(Pop_Factor), pch = pchs)
nMDS_spp
##
## Call:
## metaMDS(comm = TransSpp, distance = "bray", k = 2)
##
## global Multidimensional Scaling using monoMDS
##
## Data: TransSpp
## Distance: bray
##
## Dimensions: 2
## Stress: 0.0969794
## Stress type 1, weak ties
## No convergent solutions - best solution after 20 tries
## Scaling: centring, PC rotation, halfchange scaling
## Species: expanded scores based on 'TransSpp'
plot_spp1 <- ordiplot(nMDS_spp, display = "species", type = "n")
points(nMDS_spp, col = "black", pch = pchs[Pop_Factor])
legend("topright", bty = "n", legend = levels(Pop_Factor), pch = pchs)
nMDS_spp
##
## Call:
## metaMDS(comm = TransSpp, distance = "bray", k = 2)
##
## global Multidimensional Scaling using monoMDS
##
## Data: TransSpp
## Distance: bray
##
## Dimensions: 2
## Stress: 0.0969794
## Stress type 1, weak ties
## No convergent solutions - best solution after 20 tries
## Scaling: centring, PC rotation, halfchange scaling
## Species: expanded scores based on 'TransSpp'
plotinc <- ordiplot(spp.mds0, type = "t")
## species scores not available
pchs <- c(0:5)
Income_Factor <- factors$Income
nMDS_spp<- metaMDS(TransSpp, distance="bray", k=2)
## Run 0 stress 0.0969794
## Run 1 stress 0.1266484
## Run 2 stress 0.1266776
## Run 3 stress 0.1236921
## Run 4 stress 0.1270373
## Run 5 stress 0.09910931
## Run 6 stress 0.1236998
## Run 7 stress 0.134555
## Run 8 stress 0.1066722
## Run 9 stress 0.1198341
## Run 10 stress 0.1266631
## Run 11 stress 0.09911034
## Run 12 stress 0.1284552
## Run 13 stress 0.1051754
## Run 14 stress 0.1260333
## Run 15 stress 0.09978155
## Run 16 stress 0.1076054
## Run 17 stress 0.1075795
## Run 18 stress 0.1062512
## Run 19 stress 0.1055391
## Run 20 stress 0.09697946
## ... Procrustes: rmse 3.860572e-05 max resid 0.0001006442
## ... Similar to previous best
## *** Solution reached
plot_sppincome <- ordiplot(nMDS_spp, display = "species", type = "n")
points(nMDS_spp, col = "black", pch = pchs[Income_Factor])
legend("topright", bty = "n", legend = levels(Income_Factor), pch = pchs)
adonis(TransSpp ~ Income*Population, data = factors, permutations = 999)
##
## Call:
## adonis(formula = TransSpp ~ Income * Population, data = factors, permutations = 999)
##
## Permutation: free
## Number of permutations: 999
##
## Terms added sequentially (first to last)
##
## Df SumsOfSqs MeanSqs F.Model R2 Pr(>F)
## Income 1 0.02536 0.025357 0.50286 0.01492 0.697
## Population 1 0.02312 0.023118 0.45845 0.01360 0.744
## Income:Population 1 0.03767 0.037667 0.74697 0.02216 0.516
## Residuals 32 1.61366 0.050427 0.94932
## Total 35 1.69980 1.00000
No significant difference between groups.
library(readr)
simp <- read_csv("simp.csv", col_types = cols(Population = col_factor(levels = c("LOW",
"HIGH")), Income = col_factor(levels = c("LOW",
"HIGH"))))
```r
shapiro.test(simp$SimpsonsD)
##
## Shapiro-Wilk normality test
##
## data: simp$SimpsonsD
## W = 0.90166, p-value = 0.003789
shapiro.test(simp$Eveness)
##
## Shapiro-Wilk normality test
##
## data: simp$Eveness
## W = 0.96847, p-value = 0.3852
simp.aov <- aov(SimpsonsD ~ Income*Population, data = simp)
summary(simp.aov)
## Df Sum Sq Mean Sq F value Pr(>F)
## Income 1 0.00208 0.002079 0.805 0.376
## Population 1 0.00052 0.000520 0.201 0.657
## Income:Population 1 0.00741 0.007413 2.869 0.100
## Residuals 32 0.08270 0.002584
No significant interaction or main effect on diversity.