forked from geocompx/geocompr
-
Notifications
You must be signed in to change notification settings - Fork 0
/
_11-ex.Rmd
103 lines (82 loc) · 5.83 KB
/
_11-ex.Rmd
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
```{asis 11-ex-asis1, message=FALSE}
The solutions assume the following packages are attached (other packages will be attached when needed):
```
```{r setup11, include=FALSE}
poly_centroid = function(poly_mat) {
Origin = poly_mat[1, ] # create a point representing the origin
i = 2:(nrow(poly_mat) - 2)
T_all = lapply(i, function(x) {rbind(Origin, poly_mat[x:(x + 1), ], Origin)})
C_list = lapply(T_all, t_centroid)
C = do.call(rbind, C_list)
A = vapply(T_all, t_area, FUN.VALUE = double(1))
c(weighted.mean(C[, 1], A), weighted.mean(C[, 2], A))
}
```
```{r 11-ex-e0}
library(sf)
```
E1. Read the script [`11-centroid-alg.R`](https://github.com/geocompx/geocompr/blob/main/code/11-centroid-alg.R) in the `code` folder of the book's GitHub repository.
- Which of the best practices covered in Section \@ref(scripts) does it follow?
- Create a version of the script on your computer in an IDE\index{IDE} such as RStudio\index{RStudio} (preferably by typing the script line-by-line, in your own coding style and with your own comments, rather than copy-pasting --- this will help you learn how to type scripts). Using the example of a square polygon (e.g., created with `poly_mat = cbind(x = c(0, 9, 9, 0, 0), y = c(0, 0, 9, 9, 0))`) execute the script line-by-line.
- What changes could be made to the script to make it more reproducible?
- How could the documentation be improved?
```{asis 11-ex-e1, message=FALSE}
The script is stored in a logical location with a sensible file name.
The script is well documented with comments and the code is well formatted.
The script is reproducible.
Open a file and create a new script in RStudio, e.g., with the keyboard shortcut `Ctrl + Shift + N` (Windows) or `Cmd + Shift + N` (Mac), by clicking `File > New File > R Script` or by clicking the `+` icon in the top left of the `Source` pane.
You can also create a new R script from the R console with the command `file.create("11-centroid-alg.R")`.
The script is already reproducible, with a message stating that it needs an object called `poly_mat` to be present and, if none is present, it creates an example dataset at the outset for testing.
For people new to R it could also contain a comment stating that R must be installed before running the script.
Documentation could be improved with a more detailed description of the algorithm, including a link to the relevant section of the book.
Furthermore, the anonymous functions could be replaced with named functions and documented with Roxygen2 comments.
```
E2. In the geometric algorithms section, we calculated that the area of the polygon `poly_mat` was 245 units squared and that its centroid as at the coordinates (8.8, 9.2).
- Reproduce the results on your own computer with reference to the script [`11-centroid-alg.R`](https://github.com/geocompx/geocompr/blob/main/code/11-centroid-alg.R), an implementation of this algorithm (bonus: type the commands - try to avoid copy-pasting).
- Are the results correct? Verify them by converting `poly_mat` into an `sfc` object (named `poly_sfc`) with `st_polygon()` (hint: this function takes objects of class `list()`) and then using `st_area()` and `st_centroid()`.
```{r 11-ex-e}
# We can verify the answer by converting `poly_mat` into a simple feature collection
# as follows, which shows the calculations match:
x_coords = c(10, 20, 12, 0, 0, 10)
y_coords = c(0, 15, 20, 10, 0, 0)
poly_mat = cbind(x_coords, y_coords)
poly_sfc = sf::st_polygon(list(poly_mat))
sf::st_area(poly_sfc)
sf::st_centroid(poly_sfc)
# By calling the script:
# source("https://github.com/geocompx/geocompr/raw/main/code/11-centroid-alg.R")
```
E3. It was stated that the algorithm\index{algorithm} we created only works for *convex hulls*. Define convex hulls\index{convex hull} (see the geometry operations chapter) and test the algorithm on a polygon that is *not* a convex hull.
```{r 11-ex-e3}
x_coords = c(10, 20, 12, 0, 0, 5, 10)
y_coords = c(0, 15, 20, 10, 0, 5, 0)
plot(x_coords, y_coords, type = "l")
poly_mat = cbind(x_coords, y_coords)
# source("https://github.com/geocompx/geocompr/raw/main/code/11-centroid-alg.R")
# Area from our script: 270
poly_sfc = sf::st_polygon(list(poly_mat))
sf::st_area(poly_sfc) # Actual area: 220
```
- Bonus 1: Think about why the method only works for convex hulls and note changes that would need to be made to the algorithm to make it work for other types of polygon.
- Bonus 2: Building on the contents of `11-centroid-alg.R`, write an algorithm only using base R functions that can find the total length of linestrings represented in matrix form.
<!-- Todo: add example of matrix representing a linestring, demonstrate code to verify the answer, suggest alternative functions to decompose as a bonus. -->
```{asis 11-ex-e3-bonus1}
The algorithm would need to be able to have negative as well as positive area values.
We leave Bonus 2 as an exercise for the reader.
```
E4. In the functions section, we created different versions of the `poly_centroid()` function that generated outputs of class `sfg` (`poly_centroid_sfg()`) and type-stable `matrix` outputs (`poly_centroid_type_stable()`).
Further extend the function by creating a version (e.g., called `poly_centroid_sf()`) that is type stable (only accepts inputs of class `sf`) *and* returns `sf` objects (hint: you may need to convert the object `x` into a matrix with the command `sf::st_coordinates(x)`).
- Verify if it works by running `poly_centroid_sf(sf::st_sf(sf::st_sfc(poly_sfc)))`
- What error message do you get when you try to run `poly_centroid_sf(poly_mat)`?
```{r 11-ex-e4}
poly_centroid_sf = function(x) {
stopifnot(is(x, "sf"))
xcoords = sf::st_coordinates(x)
centroid_coords = poly_centroid(xcoords)
centroid_sf = sf::st_sf(geometry = sf::st_sfc(sf::st_point(centroid_coords)))
centroid_sf
}
poly_centroid_sf(sf::st_sf(sf::st_sfc(poly_sfc)))
poly_centroid_sf(poly_sfc)
poly_centroid_sf(poly_mat)
```