-
Notifications
You must be signed in to change notification settings - Fork 0
/
Copy pathREADME.Rmd
788 lines (661 loc) · 61.3 KB
/
README.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
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245
246
247
248
249
250
251
252
253
254
255
256
257
258
259
260
261
262
263
264
265
266
267
268
269
270
271
272
273
274
275
276
277
278
279
280
281
282
283
284
285
286
287
288
289
290
291
292
293
294
295
296
297
298
299
300
301
302
303
304
305
306
307
308
309
310
311
312
313
314
315
316
317
318
319
320
321
322
323
324
325
326
327
328
329
330
331
332
333
334
335
336
337
338
339
340
341
342
343
344
345
346
347
348
349
350
351
352
353
354
355
356
357
358
359
360
361
362
363
364
365
366
367
368
369
370
371
372
373
374
375
376
377
378
379
380
381
382
383
384
385
386
387
388
389
390
391
392
393
394
395
396
397
398
399
400
401
402
403
404
405
406
407
408
409
410
411
412
413
414
415
416
417
418
419
420
421
422
423
424
425
426
427
428
429
430
431
432
433
434
435
436
437
438
439
440
441
442
443
444
445
446
447
448
449
450
451
452
453
454
455
456
457
458
459
460
461
462
463
464
465
466
467
468
469
470
471
472
473
474
475
476
477
478
479
480
481
482
483
484
485
486
487
488
489
490
491
492
493
494
495
496
497
498
499
500
501
502
503
504
505
506
507
508
509
510
511
512
513
514
515
516
517
518
519
520
521
522
523
524
525
526
527
528
529
530
531
532
533
534
535
536
537
538
539
540
541
542
543
544
545
546
547
548
549
550
551
552
553
554
555
556
557
558
559
560
561
562
563
564
565
566
567
568
569
570
571
572
573
574
575
576
577
578
579
580
581
582
583
584
585
586
587
588
589
590
591
592
593
594
595
596
597
598
599
600
601
602
603
604
605
606
607
608
609
610
611
612
613
614
615
616
617
618
619
620
621
622
623
624
625
626
627
628
629
630
631
632
633
634
635
636
637
638
639
640
641
642
643
644
645
646
647
648
649
650
651
652
653
654
655
656
657
658
659
660
661
662
663
664
665
666
667
668
669
670
671
672
673
674
675
676
677
678
679
680
681
682
683
684
685
686
687
688
689
690
691
692
693
694
695
696
697
698
699
700
701
702
703
704
705
706
707
708
709
710
711
712
713
714
715
716
717
718
719
720
721
722
723
724
725
726
727
728
729
730
731
732
733
734
735
736
737
738
739
740
741
742
743
744
745
746
747
748
749
750
751
752
753
754
755
756
757
758
759
760
761
762
763
764
765
766
767
768
769
770
771
772
773
774
775
776
777
778
779
780
781
782
783
784
785
786
787
788
---
output: github_document
---
<!-- README.md is generated from README.Rmd. Please edit that file -->
```{r, include = FALSE}
knitr::opts_chunk$set(
collapse = TRUE,
comment = "#>",
fig.path = "man/figures/README-",
out.width = "100%"
)
```
# `plantTracker` <a href="http://www.astearsresearch.com/package/planttracker/"><img src="man/figures/logo.png" align="right" height="139" /></a>
<!-- badges: start -->
[![CRAN status](https://www.r-pkg.org/badges/version/plantTracker)](https://CRAN.R-project.org/package=plantTracker)
[![check-standard](https://github.com/aestears/plantTracker/actions/workflows/check-standard.yaml/badge.svg)](https://github.com/aestears/plantTracker/actions/workflows/check-standard.yaml)
[![test-coverage](https://github.com/aestears/plantTracker/actions/workflows/test-coverage.yaml/badge.svg)](https://github.com/aestears/plantTracker/actions/workflows/test-coverage.yaml)
[![License: MIT](https://img.shields.io/badge/License-MIT-yellow.svg)](https://opensource.org/license/mit/)
[![DOI](https://zenodo.org/badge/353094728.svg)](https://zenodo.org/badge/latestdoi/353094728)
<!-- badges: end -->
Welcome to `plantTracker`! This package was designed to transform **long-term quadrat maps** that show plant occurrence and size into **demographic data** that can be used to answer questions about population and community ecology.
## Table of Contents
* [Installing `plantTracker`](#Installation)
* [Contributing](#issues)
* [License](#license)
* [Contact](#contact)
* [Interactive tools to explore plantTracker functions](#ShinyApps)
* [How to Use `plantTracker`](#how_to_use)
* [Prepare data for use in `plantTracker`](#prep_data)
* [Ensure data is in the correct format with `checkDat()` ](#check_dat)
* [Track individuals through time using `trackSpp()`](#trackSpp)
* [Calculate local neighborhood density using `getNeighbors()`](#getNeigh)
* [Further analysis with more `plantTracker` functions](#next)
<a id="Installation"></a>
## Installation
Install `plantTracker` from
[CRAN](https://CRAN.R-project.org/package=plantTracker):
``` r
install.packages("plantTracker")
```
Alternatively, you can install the current version of `plantTracker` from
[GitHub](https://github.com/):
``` r
install.packages("devtools")
devtools::install_github("aestears/plantTracker")
```
<a id="issues"></a>
## Contributing
Please report any problems that you encounter while using `plantTracker` as "issues" on (our GitHub repository)[https://github.com/aestears/plantTracker/issues/]. Help us make this package better!
<a id="license"></a>
## License
This package is licensed under MIT License
Copyright (c) 2022 Alice Stears
<a id="contact"></a>
## Contact
Questions about `plantTracker` can be forwarded to Alice Stears, the package
maintainer, at [email protected].
<a id="ShinyApps"></a>
## Interactive tools to explore plantTracker functions
Follow this link to a shiny app that lets you change arguments in the trackSpp() function, and see how that impacts individual ID assignment in an example dataset: [https://astears.shinyapps.io/shinyapp/](https://astears.shinyapps.io/shinyapp/)
<a id="how_to_use"></a>
## How to use the `plantTracker` R package
The material below explains how to use `plantTracker`, starting with formatting your data correctly. This information is also available in the 'Suggested `plantTracker` Workflow' vignette, which is included in the package.
```{r setup, include = FALSE}
library(plantTracker)
library(ggplot2)
library(sf)
```
<a id="prep_data"></a>
### *1.* Prepare data
The functions in `plantTracker` require data in a specific format. `plantTracker` includes an example dataset that consists of two pieces: `grasslandData` and `grasslandInventory`. You can load these example datasets into your global environment by calling `data(grasslandData)` and `data(grasslandInventory)`. You can view the documentation for these datasets by calling`?grasslandData` and `?grasslandInventory`.
Most `plantTracker` functions require two data objects. The first is a data frame that contains the location and metadata for each mapped individual, which we from now on will call `dat`. The second is a list that contains a vector of years in which each quadrat was sampled, which we from now on will cal `inv`.
Below are the basic requirements for these data objects.
<a id="dat_data"></a>
#### *1.1* The `dat` data frame must . . .
* ... be an `sf` data.frame. More on this below in section [*1.1.1*](#dat_to_sf)...
* ... contain a row for each individual observation in each year.
* ... have a column that contains character strings indicating the specific epithet for each observation. This column must be a character data type. The function expects this column to be called "Species", but a different name can be specified in function calls.
* ... have a column that contains character strings indicating the site at which each observation was collected. This is a level of classification "above" the quadrat (i.e. all quadrats measured at the Central Plains Experimental Range in Nunn, CO might have the value "CO" in this column). This column must be a character data type. The function expects this column to be called "Site", but a different name can be specified in function calls.
* ... have a column that contains character strings indicating the quadrat at which each observation was collected. This column must be a character data type. The function expects this column to be called "Quad", but a different name can be specified in function calls.
* ... have a column that contains a value indicating the year when this individual observation was collected. This must be a numeric data type, and must be either a four or two digit year number. The function expects this column to be called "Year", but a different name can be specified in function calls.
* ... have a column (almost always called "geometry" in the `sf` package data format) that contains a polygon representing the location of each observation. Each observation must be a `POLYGON` or `MULTIPOLYGON`. Data cannot be stored as `POINTS`.
+ If the data was collected such that forbs or small grasses were mapped as points and digitized as such, then these observations must be converted to polygons. We recommend that you convert them to small circular polygons with an identical radius. If you do this transformation, we also recommend that you include a column that indicates whether each row was originally mapped as a polygon or a point, since the demographic data that deals with size will be relatively meaningless for observations originally mapped as points.
+ `dat` does not need to have a coordinate reference system (i.e. CRS can be "NA"), but it can have one if you'd like.
* ... *not* have columns called "neighbors", "nearEdge", "trackID", "age", "size_tplus1", "recruit", "survives_tplus1", "basalArea_ramet", or "basalArea_genet", since these columns are added by `plantTracker` functions.
* Note: Additional columns can be included in the input data.frame, although
they might not be included in the output of `plantTracker` functions.
Here are the first few rows of a possible `dat` input data.frame:
```{r echo = FALSE}
dataShow <- head(grasslandData[, !names(grasslandData) %in%
c("Clone", "Seedling", "Stems", "Basal",
"sp_code_4", "Area")])
rownames(dataShow) <- NULL
(dataShow)
#knitr::kable(dataShow, caption = "**Table 1.1**: Example `dat` data.frame" )
```
* Note: that the required columns are "Species", "Site", "Quad", "Year", and "geometry". The additional columns "Type" and "sp_code_6" are just "along for the ride" in any analysis using `plantTracker` functions.
Here's what some of the example `dat` data (from the "SG2" quadrat at the "AZ" site in 1922) look like when plotted spatially:
```{r, echo = FALSE, fig.width=7, fig.align = 'center', fig.cap = "Figure 1.1 : Spatial map of a subset of example 'dat' dataset"}
exampleDat <- grasslandData[grasslandData$Site == "AZ" &
grasslandData$Quad == "SG2" &
grasslandData$Year == 1922, ]
dat <- grasslandData[grasslandData$Site == "AZ" &
grasslandData$Quad == "SG2",]
ggplot(data = exampleDat) +
geom_sf(aes(color = Species, fill = Species)) +
geom_segment(aes(x = 0, xend = 1, y = 0, yend = 0), size = .5,
lineend = "round", color = "grey30") +
geom_segment(aes(x = 0, xend = 1, y = 1, yend = 1), size = .5,
lineend = "round", color = "grey30") +
geom_segment(aes(x = 0, xend = 0, y = 0, yend = 1), size = .5,
lineend = "round", color = "grey30") +
geom_segment(aes(x = 1, xend = 1, y = 0, yend = 1), size = .5,
lineend = "round", color = "grey30") +
xlab("quadrat horizontal edge (m)") +
ylab("quadrat vertical edge (m)") +
theme_classic() +
theme(axis.line = element_blank(),
legend.text = element_text(face = "italic"),
plot.margin = margin(1,0,1,0))
```
It's important to note that, while `plantTracker` was designed to be used with small-scale maps of plant occurrence in quadrats, it is conceivably possible to use other styles of map data in `plantTracker` functions. All that is required is a single mapped basal area (or point location converted to a small polygon) at each time point for each organism (or ramet), and is accompanied by the required metadata detailed above. For example, `plantTracker` functions could be used to estimate tree demographic rates at the scale of 100 m x 50 m plots.
<a id="dat_to_sf"></a>
##### *1.1.1* Get your data into the `sf` data format
As mentioned above, `plantTracker` uses the `sf` R package to deal with spatial
data. The map data that `plantTracker` was built to analyze is inherently spatial,
so you need to know how to the basics of dealing with spatial data in R if you
want to use `plantTracker`! There are many good resources to help you orient
yourself to working with spatial data in R generally:
* [https://cengel.github.io/R-spatial/intro.html](https://cengel.github.io/R-spatial/intro.html)
* [https://www.r-bloggers.com/2021/06/using-geospatial-data-in-r/](https://www.r-bloggers.com/2021/06/using-geospatial-data-in-r/)
And the `sf` package more specifically:
* [https://CRAN.R-project.org/package=sf](https://CRAN.R-project.org/package=sf)
* [https://r-spatial.github.io/sf/](https://r-spatial.github.io/sf/)
These resources provide a great orientation, and while I recommend looking over
them if you're new to working with spatial data in R, I've included a brief
tutorial for uploading shapefiles into R as `sf` data frames.
Most of the published chart-quadrat datasets have the map data stored as
shapefiles in complex file structures, which can be a bit confusing to navigate.
`plantTracker` requires all of your data (for all species, plots and years) to be
in one single data frame. *This example shows how you might navigate through a
complex file structure to to pull out shapefiles and put them into one single
`sf` data frame. for further analysis with `plantTracker`.*
For this example, I'll use a subset of the data from the Santa Rita Experimental
Range in Arizona, which has been published in [this data paper](https://figshare.com/articles/dataset/Data_Paper_Data_Paper/3553581?file=5622201). In this dataset, shapefiles for each quadrat are stored in their own folder.
Within that folder there are two shapefiles for each year: one that contains map
data for polygons, and one that contains data for points. The following code
reads in those shapefiles, transforms the points to polygons of a fixed radius,
and puts all the data into one `sf` data frame. If you want to follow along,
download the "shapefiles.zip" file from the data paper, un-zip it, and name it "AZ_shapefiles". The dataset that is the result of this example is the same as
part of the "grasslandData" dataset included in `plantTracker`.
```{r, echo=FALSE, include=FALSE, eval = FALSE}
# get the working directory name
wdName <- c("/Users/Alice/Dropbox/Grad School/Research/Trait Project/Data/ChartQuadDatasets/AZs_downloaded data/")
```
```{r, warning=FALSE, eval = FALSE}
# save a character vector of the file names in the file that contains the
# shapefiles (in this case, called "CO_shapefiles"), each of which is a quadrat
# note: 'wdName' is a character string indicating the path of the directory
# containing the 'AZ_shapefiles' folder
quadNames <- list.files(paste0(wdName,"AZ_shapefiles/"))
# trim the quadrats down to 2, for the sake of runtime in this example
quadNames <- quadNames[quadNames %in% c("SG2", "SG4")]
# now we'll loop through the quadrat folders to download the data
for (i in 1:2){#length(quadNames)) {
# get the names of the quadrat for this iteration of the loop
quadNow <- quadNames[i]
# get a character vector of the unique quad/Year combinations of data in
# this folder that contain polygon data
quadYears <- quadYears <- unlist(strsplit(list.files(
paste0(wdName, "AZ_shapefiles/",quadNow,"/"),
pattern = ".shp$"), split = ".shp"))
# loop through each of the years in this quadrat
for (j in 1:length(quadYears)) {
# save the name of this quadYear combo
quadYearNow <- quadYears[j]
# read in the shapefile for this quad/year combo as an sf data frame
# using the 'st_read()' function from the sf package
shapeNow <- sf::st_read(dsn = paste0(wdName,"AZ_shapefiles/",quadNow),
# the 'dsn' argument is the folder that
# contains the shapefile files--in this case,
# the folder for this quadrat
layer = quadYearNow) # the 'layer' argument has the
# name of the shapefile, without the filetype extension! This is because each
# shapefile consists of at least three separate files, each of which has a
# unique filetype extension.
# the shapefiles in this dataset do not have all of the metadata we
# need, and have some we don't need, so we'll remove what we don't need and
# add columns for 'site', 'quad', and 'year'
shapeNow$Site <- "AZs"
shapeNow$Quad <- quadNow
# get the Year for the shapefile name--in this case it is the last for
# numbers of the name
shapeNow$Year <- as.numeric(strsplit(quadYearNow, split = "_")[[1]][2]) + 1900
# determine if the 'current' quad/year contains point data or polygon data
if (grepl(quadYearNow, pattern = "C")) { # if quadYearNow has point data
# remove the columns we don't need
shapeNow <- shapeNow[,!(names(shapeNow)
%in% c("Clone", "Seedling", "Area", "Length", "X", "Y"))]
# reformat the point into a a very small polygon
# (a circle w/ a radius of .003 m)
shapeNow <- sf::st_buffer(x = shapeNow, dist = .003)
# add a column indicating that this observation was originally
# mapped as a point
shapeNow$type <- "point"
} else { # if quadYearNow has polygon data
# remove the columns we don't need
shapeNow <- shapeNow[,!(names(shapeNow) %in% c("Seedling", "Canopy_cov", "X", "Y", "area"))]
# add a column indicating that this observation was originally
# mapped as a polygon
shapeNow$type <- "polygon"
}
# now we'll save this sf data frame
if (i == 1 & j == 1) { # if this is the first year in the first quadrat
dat <- shapeNow
} else { # if this isn't the first year in the first quadrat, simply rbind
# the shapeNow sf data frame onto the previous data
dat <- rbind(dat, shapeNow)
}
}
}
# Now, all of the spatial data are in one sf data frame!
# for the sake of this example, we'll remove data for some species and years in order to make the example run faster (and to make this 'dat' data.frame identical to the "grasslandData" dataset included in this R pakcage).
dat <- dat[dat$Species %in% c("Heteropogon contortus", "Bouteloua rothrockii", "Ambrosia artemisiifolia", "Calliandra eriophylla", "Bouteloua gracilis", "Hesperostipa comata", "Sphaeralcea coccinea", "Allium textile"),]
dat <- dat[ (dat$Quad %in% c("SG2", "SG4") &
dat$Year %in% c(1922:1927)),]
```
In some spatial datasets, observations that were measured as "points" in the field are still stored as "points" in the shapefiles. `plantTracker` requires all observations to be stored as "polygon" geometry in order to streamline functions, so we need to translate "points" into small polygons of a fixed area. In this case, we'll transform them into circles with a radius of 1 cm (.01, since this dataset measures area in meters). 'dat' has a column called "type." A value of "point" in this column will tell us that, even though the geometry of the "point" data is now in "polygon" format, the values for basal area and growth are not indicative of the true size of the plant.
```{r}
# We use the function "st_buffer()" to add a buffer of our chosen radius (.01) around each point observation, which will transform each observation into a circle of the "polygon" format with a radius of .01.
dat_1 <- st_buffer(x = dat[st_is(x = dat, type = "POINT"),], dist = .01)
dat_2 <- dat[!st_is(x = dat, type = "POINT"),]
dat <- rbind(dat_1, dat_2)
```
If you don't want to download the data and format it into an sf data.frame, you
can also use a subset of the "grasslandData" data object stored in this R
package. You will just need to subset it to include only the data from the "AZ"
site. Code to do this is below:
```{r}
dat <- grasslandData[grasslandData$Site == "AZ",]
```
<a id="dat_inv"></a>
#### *1.2* The `inv` list must . . .
* ... be a named list
* ... have element names that are each a character string identical to the name of a quadrat in `dat`. There cannot be two elements with the same name, and there cannot be an element with more than one quadrat in its name. There must be an element for each quadrat with data in `dat`.
* ... have list elements that are each a numeric vector of years in which each quadrat was sampled. If a quadrat was "skipped" in a year, then that year must be excluded from this vector. The format of the years must be the same as the format of years in `dat` (i.e. if year is a four-digit number in `dat`, then it must be a four-digit number in `inv`). Make sure this is the years the quadrat was actually sampled, not just the years that have data in the `dat` data frame! This argument allows the function to differentiate between years when the quadrat wasn't sampled and years when there just weren't any individuals of a species present in that quadrat. If a quadrat wasn't sampled in a given year, don't put an 'NA' in `inv` for that year! Instead, just skip that year.
Here is an example of an `inv` argument that corresponds to the example `dat` argument above. The quadrats that have data in `dat` are "SG2" and "SG4", so there are elements in `inv` that correspond to each of these quadrats.
```{r echo=FALSE, fig.align='center', fig.cap = "Figure 1.2: An example of the `inv` argument"}
(inv <- grasslandInventory[c( "SG2", "SG4")])
```
If you already have a quadrat inventory as a data frame, it isn't complicated to reformat it to work with `plantTracker` functions.
For example, if your quadrat inventory data frame looks like this... :
```{r, echo = FALSE}
data.frame("quad1" = c(2000, 2001, NA, 2003, 2004, 2005, 2006, 2007),
"quad2" = c(2000:2007),
"quad3" = c(2000, NA, 2002, 2003, 2004, 2005, 2006, 2007))
```
... then do the following to get it into a format ready for `plantTracker`:
```{r}
quadInv_DF <- data.frame("quad1" = c(2000, 2001, NA, 2003, 2004, 2005, 2006, 2007),
"quad2" = c(2000:2007),
"quad3" = c(2000, NA, 2002, 2003, 2004, 2005, 2006, 2007))
# use the 'as.list()' function to transform your data frame into a named list
quadInv_list <- as.list(quadInv_DF)
# we still need to remove the 'NA' values, which we can do using the
# 'lapply()' function
(quadInv_list <- lapply(X = quadInv_list, FUN = function(x) x[is.na(x) == FALSE]))
```
<a id="check_dat"></a>
#### *1.3* Check the `inv` and `dat` arguments using `checkDat()`
The generic `checkDat()` function:
```r
checkDat(dat, inv = NULL, species = "Species", site = "Site", quad = "Quad",
year = "Year", geometry = "geometry", reformatDat = FALSE, ...)
```
This step is optional, but can be useful if you're unsure whether your `dat` and `inv` arguments are in the correct format. The `plantTracker` function `checkDat()` takes `dat` and `inv` as arguments for the arguments `dat` and`inv`, and will return informative error messages if either argument is not in the correct format.
Additional optional arguments to `checkDat()` are `species`, `site`, `quad`, `year`, `geometry`, and `reformatDat`.
* **`species/site/quad/year/geometry`** These arguments only need to be included if the columns in `dat` that contain the data for species, site, quadrat, year and geometry of each observation are *different* from the names "Species", "Site", "Quad", "Year, and "geometry". For example, if the column in your version of `dat` that contains the species identity of each observation is called "species_names", then the argument `species = "species_names"` must be included in your call to `checkDat()`.
* **`reformatDat`** is a TRUE/FALSE argument that determines whether you want the `checkDat()` function to return a version of `dat` that is ready to go into the steps of this workflow. If `reformatDat` = TRUE then `checkDat()` will return a list that contains the reformatted version of `dat`, the reformatted version of `inv` and an additional element called "userColNames", which contains the column names in the input version of `dat` that are different from the expected column names of "Species", "Site", "Quad", "Year, and "geometry" (if there are any). If `reformatDat` = TRUE, then `checkDat()` will return a message indicating that
your data is ready for the next step. The default value is FALSE.
___
<a id="trackSpp"></a>
### *2* Track individuals through time using `trackSpp()`
Now it's time to transform your raw dataset into demographic data! This is accomplished using the `trackSpp()` function. This function follows individual plants from year to year in the same quadrat to determine survival, size in the next year, age, and some additional potentially-useful demographic data. It does this by comparing quadrat maps from sequential years. If there is overlap of individuals of the same species in consecutive years, then the rows in `dat` that contain data for those overlapping individuals are given the same "trackID", or unique identifier.
Here is the generic `trackSpp()` function:
```r
trackSpp(dat, inv, dorm, buff, buffGenet, clonal, species = "Species",
site = "Site", quad = "Quad", year = "Year", geometry = "geometry",
aggByGenet = TRUE, printMessages = TRUE, flagSuspects = FALSE,
shrink = 0.1, dormSize = 0.05, ...)
```
<a id="trackSpp_args"></a>
#### *2.1* Function arguments
`trackSpp()` takes the following arguments:
* **`dat`** This is the `sf` data frame that we've been calling `dat` so far. This must be in the correct format (which you can check before-hand using `checkDat()`), but informative error messages will be returned if it is incorrect. It *must* have the columns outlined in Section [1.1](#dat_to_sf), but they can have different names as long as those names are included in this function call (more on that later...).
* **`inv`** This is the list of quadrat sampling years we've been calling `inv`. If it is not in the correct format or does not contain data for the correct quadrats, then an informative error message will be returned.
* **`dorm`** This is a positive integer value that indicates how long you want the function to allow an individual to be "dormant". In this case, dormancy can be interpreted as the biological phenomenon where a plant has above-ground tissue present in year 1, is alive underground but with no above-ground tissue in year 2, and then has above-ground tissue in a subsequent year. Dormancy can also be interpreted here as data-collection error, whereby an individual is accidentally not mapped in between years where it was recorded.
*Consider the following example*: There is a polygon of species "A" in year 1, which is our "focal individual". In year 2, there is not a polygon of species "A" that overlaps with our focal individual. In year 3, there is a polygon of species "A" that is in the same location as our focal individual. If `dorm = 0`, then our focal individual would get a 0 in the survival column, and the polygon of species "A" in year 3 would be considered a new recruit and get a new trackID. If `dorm = 1`, because there is overlap between two polygons of the same species with only a 1-year gap between when they occur, these two polygons will be considered the same genetic individual, will have the same trackID, and our focal individual will have a "1" in the survival column. In an alternative scenario, in years 3 and 4 there are not polygons of species "A" that are in the same location as our focal individual, but there is a polygon in year 4 that overlaps our focal individual. If `dorm = 1`, then our focal individual would get a "0" for survival, but if `dorm = 2`, then it would get a 1 for survival.
```{r echo = FALSE, fig.width=8, fig.height=3, fig.align = 'center', fig.cap = "Figure 2.1: A visualization of the 'dormancy' scenario described above."}
exampleSmall <- exampleDat[unique(c(30, 42, 44, 25, 2, 61, 59, 52, 45, 8, 42, 37, 60, 34, 45, 38, 26, 49, 52, 46)),]
exampleSmall$FocalInd <- 'other individual'
exampleSmall[10,'FocalInd'] <- 'focal individual'
exampleSmall$Year <- "Year 1"
exampleSmall_2 <- exampleSmall[c(1:9,11:nrow(exampleSmall)),]
exampleSmall_2$Year <- "Year 2"
exampleSmall_3 <- exampleSmall
exampleSmall_3$Year <- "Year 3"
exampleSmall <- rbind(exampleSmall, exampleSmall_2, exampleSmall_3)
ggplot(data = exampleSmall) +
geom_sf(aes(color = FocalInd, fill = FocalInd)) +
geom_segment(aes(x = 0, xend = .6, y = 0, yend = 0), size = .5,
lineend = "round", color = "grey30") +
geom_segment(aes(x = 0, xend = .6, y = .6, yend = .6), size = .5,
lineend = "round", color = "grey30") +
geom_segment(aes(x = 0, xend = 0, y = 0, yend = .6), size = .5,
lineend = "round", color = "grey30") +
geom_segment(aes(x = .6, xend = .6, y = 0, yend = .6), size = .5,
lineend = "round", color = "grey30") +
xlab("quadrat horizontal edge (m)") +
ylab("quadrat vertical edge (m)") +
xlim(c(0,.6)) +
ylim(c(0,.6)) +
#labs(title = Year) +
facet_wrap( ~ Year) +
theme_classic() +
theme(axis.line = element_blank(),
legend.text = element_text(face = "italic"),
plot.margin = margin(.5,0,.5,0),
legend.title = element_blank()
) +
scale_fill_discrete(type = c("#E69F00", "#A6A6A6", "#009E73")) +
scale_color_discrete(type = c("#E69F00", "#A6A6A6", "#009E73"))
```
If you'd like to be more specific and perhaps biologically accurate, you can also specify the `dorm` argument uniquely for each species. For example, it might be that you are confident that your data collectors did not accidentally "miss" any individuals, and your `dat` data frame contains observations for shrubs or trees, which are very unlikely to go dormant, and small forbs, which are much more likely to go dormant for one or two years. In order to disallow dormancy for trees and shrubs, but to allow dormancy for forbs, you will provide a data frame to the `dorm` argument instead of a single positive integer value. There will be two columns: 1) a "Species" column that has the species name for each species present in `dat`, and 2) a column called "dorm" that has positive integer values indicating the dormancy you'd like to allow for each species. Make sure that if you are following the data frame approach, you must provide a dormancy argument for *every* species that has data in `dat`. Make sure that the species names in the `dorm` data frame are spelled exactly the same as they are in `dat`. The data frame should look something like this:
```{r echo = FALSE}
dormDF <- data.frame("Species" = c("tree A", "shrub B", "tree C", "forb D", "forb E", "forb F"),
"dorm" = c(0,0,0,1,2,1))
(dormDF)
# knitr::kable(dormDF, caption = "**Table 2.1**: Example 'dorm' data.frame", )
```
+ *Important Note:* Be very careful about how you define the `dorm` argument. The bigger the `dorm` argument, the more likely you are to overestimate survival. For annually-sampled data, I would need a very biologically-compelling reason to to specify a `dorm` argument greater than 1 year.
* **`buff`** This is a positive numeric value that indicates how much an individual can move from year 1 to year 2 and still be considered the same individual (receive the same trackID). In addition to accounting for true variation in location of a plant's stem from year to year, this argument also accounts for small inconsistencies in mapping from year to year. The `buff` argument must be in the same units as the spatial values in `dat`. For example, if the spatial data in `dat` is measured in meters, and you want to allow a plant to "move" 15 cm between year 1 and year 2, then you would include the argument `buff = .15` in your call to `trackSpp()`. If you want to allow no movement, use `buff = 0`. Below is a visualization of two different `buff` scenarios.
```{r echo = FALSE, warning= FALSE, fig.width=7, fig.height = 3, fig.align = 'center', fig.cap = "Figure 2.2: With a 10 cm buffer, these polygons in 1922 and 1923 overlap and will be identified by trackSpp() as the **same** individual and receive the same trackID."}
exampleDat <- grasslandData[grasslandData$Site == "AZ" &
grasslandData$Quad == "SG2" &
grasslandData$Year %in% c(1922, 1923),]
exampleDatIDsTemp <-
exampleDat[exampleDat$Species == "Bouteloua rothrockii", ]
exampleDatIDsTemp <-
exampleDatIDsTemp[round(exampleDatIDsTemp$Area, 7) %in% round(c(0.0005471808, 0.0005321236), 7), ]
exampleDatIDsTemp$ghost <- "observation from current year"
exampleBuffed <-
st_buffer(exampleDatIDsTemp[round(exampleDatIDsTemp$Area, 7) == 0.0005472, ], dist = .10)
exampleBuffed$Year <- 1922
exampleBuffed$ghost <- "10 cm buffer"
exampleBuffedNext <- exampleBuffed
exampleBuffedNext$Year <- 1923
ghost <-
exampleDatIDsTemp[round(exampleDatIDsTemp$Area, 7) == 0.0005472, ]
ghost$Year <- 1923
ghost$ghost <- "polygon location in previous year"
exampleDatIDs <-
rbind(exampleDatIDsTemp, ghost, exampleBuffed, exampleBuffedNext)
ggplot(data = exampleDatIDs) +
geom_sf(aes(
fill = ghost,
alpha = ghost,
color = ghost,
lty = ghost
)) +
geom_segment(
aes(
x = .6,
xend = 1,
y = .4,
yend = .4
),
size = .5,
lineend = "round",
color = "grey30"
) +
geom_segment(
aes(
x = .6,
xend = 1,
y = .8,
yend = .8
),
size = .5,
lineend = "round",
color = "grey30"
) +
geom_segment(
aes(
x = .6,
xend = .6,
y = .4,
yend = .8
),
size = .5,
lineend = "round",
color = "grey30"
) +
geom_segment(
aes(
x = 1,
xend = 1,
y = .4,
yend = .8
),
size = .5,
lineend = "round",
color = "grey30"
) +
xlab("quadrat horizontal edge (m)") +
ylab("quadrat vertical edge (m)") +
#labs(title = Year) +
facet_wrap(~ Year) +
theme_classic() +
theme(
axis.line = element_blank(),
legend.text = element_text(face = "italic"),
#plot.margin = margin(1,0,1,0),
legend.title = element_blank()
) +
scale_fill_discrete(type = c("#E69F00", "#009E73", "#A6A6A6")) +
scale_color_discrete(type = c("#E69F00", "#009E73", "#A6A6A6")) +
scale_alpha_discrete(range = c(0, 1, .5)) +
scale_linetype_manual(values = c("twodash", "solid", "dotted"))
```
```{r echo = FALSE, warning= FALSE, fig.width=7, fig.height = 3, fig.align = 'center', fig.cap = "Figure 2.3: With a 3 cm buffer, these polygons in 1922 and 1923 don't quite overlap, so will be identified by trackSpp() as **different** individuals and receive different trackIDs."}
exampleBuffed <- st_buffer(exampleDatIDsTemp[round(exampleDatIDsTemp$Area, 7) ==0.0005472,], dist = .03)
exampleBuffed$Year <- 1922
exampleBuffed$ghost <- "3 cm buffer"
exampleBuffedNext <- exampleBuffed
exampleBuffedNext$Year <- 1923
ghost <- exampleDatIDsTemp[round(exampleDatIDsTemp$Area, 7) ==0.0005472,]
ghost$Year <- 1923
ghost$ghost <- "polygon location in previous year"
exampleDatIDs <- rbind(exampleDatIDsTemp, ghost, exampleBuffed, exampleBuffedNext )
ggplot(data = exampleDatIDs
) +
geom_sf(aes(fill = ghost, alpha = ghost, color = ghost, lty = ghost)) +
geom_segment(aes(x = .6, xend = 1, y = .4, yend = .4), size = .5,
lineend = "round", color = "grey30") +
geom_segment(aes(x = .6, xend = 1, y = .8, yend = .8), size = .5,
lineend = "round", color = "grey30") +
geom_segment(aes(x = .6, xend = .6, y = .4, yend = .8), size = .5,
lineend = "round", color = "grey30") +
geom_segment(aes(x = 1, xend = 1, y = .4, yend = .8), size = .5,
lineend = "round", color = "grey30") +
xlab("quadrat horizontal edge (m)") +
ylab("quadrat vertical edge (m)") +
#labs(title = Year) +
facet_wrap( ~ Year) +
theme_classic() +
theme(axis.line = element_blank(),
legend.text = element_text(face = "italic"),
#plot.margin = margin(1,0,1,0),
legend.title = element_blank()) +
scale_fill_discrete(type = c("#E69F00", "#009E73", "#A6A6A6")) +
scale_color_discrete(type = c("#E69F00", "#009E73", "#A6A6A6")) +
scale_alpha_discrete(range = c(0, 1, .5)) +
scale_linetype_manual(values=c("twodash", "solid", "dotted"))
```
* **`clonal`** This is a logical argument (TRUE or FALSE) that indicates whether you want to allow plants to be clonal or not. In the context of this type of data, "clonal" means that one genetic individual (or "genet") can be recorded as multiple polygons (or "ramets"). If `clonal = TRUE`, then multiple polygons in the same year can be part of the same individual and have the same trackID. If `clonal = FALSE`, then every polygon in a given year is a unique individual and has a unique trackID. This option can be defined globally for all species present in `dat` by setting `clonal` equal to FALSE or TRUE in the `trackSpp()` function call. Alternatively, `clonal` can be specified uniquely for each species by creating a data frame that contains a `clonal` argument for each species (analogous to the data frame for the `dorm` argument shown in Table 2.1, but with a column called "clonal").
The following arguments to `trackSpp()` are only required in certain contexts.
* **`buffGenet`** is an argument that is only required if `clonal = TRUE` or if `clonal` is a data frame that contains at least a single `TRUE` in the "clonal" column. `buffGenet` is a numeric value that indicates how close polygons of the same species must be to one another in the first year in order to be considered parts of the same genetic individual (ramets of the same genet). Similar to `buff`, be very careful and conservative when defining this argument. A large value for `buffGenet` can quickly lead to the entire quadrat being treated as the same genetic individual! I suggest experimenting with multiple values of `buffGenet`, looking at maps that show the trackID assignment, and deciding on a value that leads to trackID assignments that make the most biological sense to you. This argument is passed to the `groupByGenet()` function, which assigns the same trackID to individuals that are within the `buffGenet` buffer of each other. Polygons are only grouped by genet in the first year of data. After that point, grouping by genet happens based on data in previous years. If there are multiple polygons that overlap with a genet in the previous year, they are given the same trackID and are considered ramets belonging to the same genet. The value of `buffGenet` must be greater than or equal to zero, and must be in the same units as the spatial data in `dat`. `buffGenet` can be a single numeric value which will be applied to all species present in `dat`, or can be specified uniquely for each species by creating a data frame that contains a `buffGenet` argument for each species (analogous to the data frame for the `dorm` argument shown in Table 2.1, but with a column called "buffGenet").
* **`aggByGenet`** is a logical argument that is only required if `clonal = TRUE` or if `clonal` is a data frame that contains at least a single `TRUE` in the "clonal" column. This argument determines whether the output data frame from `trackSpp()` will have a row for every single ramet, or will be aggregated so that each genet is represented by a single row. If `aggByGenet = FALSE`, then the output is not aggregated. If `aggByGenet = TRUE` (the default setting), then the results are aggregated using the `plantTracker` function `aggregateByGenet()`. This function combines the `sf` "POLYGONS" for each ramet into one `sf` `MULTIPOLYGON` for the entire genet, and combines the associated metadata ("Species", "Site", "Quad", "Year", "trackID", "basalArea_genet", "age", "recruit", "survives_tplus1", "size_tplus1", "nearEdge") into one row for this genet. Even if the input `dat` had additional columns, they will not be included in the output of `trackSpp` if `aggByGenet = TRUE`, since it is uncertain if they can be summed across all ramets or are identical across all ramets. For example, if each ramet has a unique character string in a column called "name", there is no easy way to "sum" the character strings in this column to have one value for each genet. If you want the output data frame from `trackSpp()` to have the same columns as your input `dat` data.frame, set the `aggByGenet` argument to FALSE. However, Be Careful, since any demographic analysis should be done with a data.set that has only one row per genet, otherwise you will be estimating survival and growth rates on the scale of ramets instead of genets. If you take the `aggByGenet = FALSE` route, be sure to pass your dataset through the `aggregateByGenet()` function (or aggregate to the genet scale using your preferred method) before demographic analysis.
* **`species/site/quad/year/geometry`** These arguments only need to be included if the columns in `dat` that contain the data for species, site, quadrat, year and geometry of each observation are *different* from the names "Species", "Site", "Quad", "Year, and "geometry". For example, if the column in your version of `dat` that contains the species identity of each observation is called "species_names", then the argument `species = "species_names"` must be included in your call to `trackSpp()`.
* **`printMessages`** This is an optional logical argument that determines whether this function returns messages about genet aggregation, as well as messages indicating which year is the last year of sampling in each quadrat and which year(s) come before a gap in sampling that exceeds the `dorm` argument (and thus which years of data have an "NA" for "survives_tplus1" and "size_tplus1"). If `printMessages = TRUE` (the default), then messages are printed. If `printMessages = FALSE`, then messages are not printed.
* **`flagSuspects`** This is an optional logical argument of length 1, indicating whether observations that are "suspect" will be flagged. The default is`flagSuspects = FALSE`. If `flagSuspects = TRUE`, then a column called "Suspect" is added to the output data.frame. Any suspect observations get a "TRUE" in the "Suspect" column, while non-suspect observations receive a "FALSE". There are two ways that an observation can be classified as "suspect". First, if two consecutive observations have the same trackID, but the observation in year t+1 is less that a certain percentage (defined by the `shrink` arg.) of the observation in year t, it is possible that the observation in year t+1 is a new recruit and not the same individual. The second way an observation can be classified as "suspect" is if it is very small before going dormant. It is unlikely that a very small individual will survive dormancy, so it is possible that the function has mistakenly given a survival value of "1" to this individual. A "very small individual" is any observation with an area below a certain percentile (specified by `dormSize`) of the size distribution for this species, which is generated using all of the size data for this species in `dat`. If you are using the output dataset for demographic analysis, you may want to exclude "Suspect" observations. If `flagSuspects = FALSE`, then no additional column is added.
* **`shrink`** This is an optional argument that takes a single numeric value. This value is only used when `flagSuspects = TRUE`. When two consecutive observations have the same trackID, and the ratio of size_t+1 to size_t is smaller than the value of `shrink`, the observation in year t gets a "TRUE" in the "Suspect" column. For example, `shrink = 0.2`, and an individual that the tracking function has identified as "BOUGRA_1992_5" has an area of 9 cm$^2$ in year_t and an area of 1.35 cm$^2$ in year_t+1. The ratio of size_t+1 to size_t is 1.35/9 = 0.15, which is smaller than the cutoff specified by `shrink`, so the observation of "BOUGRA_1992_5" in year t gets a "TRUE" in the "Suspect" column. The default value is `shrink = 0.10`.
* **`dormSize`** This is an optional argument that takes a single numeric value. This value is only used when `flagSuspects = TRUE` and `dorm ≥ 1`. An individual is flagged as "suspect" if it "goes dormant" and has a size that is less than or equal to the percentile of the size distribution for this species that is designated by `dormSize`. For example, `dormSize = 0.05`, and an individual has a basal area of 0.5 cm$^2$. The 5th percentile of the distribution of size for this species, which is made using the mean and standard deviation of all observations in `dat` for the species in question, is 0.6 cm$^2$. This individual does not have any overlaps in the next year (year_t+1), but does have an overlap in year_t+2. However, because the basal area of this observation is smaller than the 5th percentile of size for this species, the observation in year t will get a "TRUE" in the "Suspect" column. It is possible that the tracking function has mistakenly assigned a "1" for survival in year_t, because it is unlikely that this individual is large enough to survive dormancy. The default value is `dormSize = .05`.
These are all of the possible arguments to `trackSpp()`!
<a id="trackSpp_out"></a>
#### *2.2* Function output
Below is an example of a potential function call to `trackSpp()`, using
the example `dat` and `inv` data we've used so far. :
```{r}
datTrackSpp <- plantTracker::trackSpp(dat = dat, inv = inv,
dorm = 1,
buff = .05,
buffGenet = .005,
clonal = data.frame("Species" = c("Heteropogon contortus",
"Bouteloua rothrockii",
"Ambrosia artemisiifolia",
"Calliandra eriophylla"),
"clonal" = c(TRUE,TRUE,FALSE,FALSE)),
aggByGenet = TRUE,
printMessages = FALSE
)
```
And here's what the output of this call to `trackSpp()` looks like:
```{r, echo = FALSE}
temp <- datTrackSpp
names(temp)[c(7,9,11)] <- c("basalArea", "survives_t+1", "size_t+1")
(temp)
#knitr::kable(head(temp[100:120,]))
```
If you did not allow any species to be clonal (`clonal = 0`) or if `aggByGenet = TRUE` in your call to `trackSpp()`, then your output data frame will have one row for each genet, and is ready for demographic analysis! If your output data frame is *not* yet aggregated by genet (i.e. you use `aggByGenet = FALSE`), then you need to transform your data frame so that each genet is represented by only one row of data. You can use the `aggregateByGenet()` function from `plantTracker` (see this function's documentation for guidance), or your own method of choice.
You can stop here and proceed to your own analyses using the demographic data you generated, or you can proceed with other `plantTracker` functions outlined below for some additional useful data.
___
<a id="getNeigh"></a>
### *3* Calculate local neighborhood density using `getNeighbors()`
It is often useful in demographic analyses to have some idea of the competition (or facilitation) that an individual organism is dealing with. Interactions between individuals can have a profound impact on whether an organism survives and grows. Spatial datasets of plant occurrence allow us to generate an estimate of the interactions an individual plant has with other plants by determining how many other individuals occupy the "local neighborhood" of each focal plant. While this isn't a direct measure of competition or facilitation, it gives us an estimate that we can include in demographic models.
Here is the generic `getNeighbors()` function:
```r
getNeighbors(dat, buff, method, compType = "allSpp", output = "summed",
trackID = "trackID", species = "Species", quad = "Quad", year = "Year",
site = "Site", geometry = "geometry", ...)
```
<a id="getNeigh_args"></a>
#### *3.1* Function options and arguments
The `getNeighbors()` function in `plantTracker` calculates local neighborhood density for each unique individual in your dataset. A user-specified buffer is drawn around each individual, and then the function counts the number of other plants within this buffer.This function can only be run on a dataset where each unique individual (genet) is represented by only one row of data. If the genet consists of multiple polygons, then they must be aggregated into one `sf` `MULTIPOLYGON` object. If your dataset has multiple rows for each genet, then you can use the `aggregateByGenet()` function to get it ready to use in `getNeighbors()`. Additionally, `getNeighbors()` requires your dataset to have a column containing a unique identifier for each genet. Across multiple years, that genet must have the same unique identifier. If you are using this function right after `trackSpp()`, your dataset will already have this unique identifier in a column called "trackID".
`getNeighbors()` has several options that allow you to customize how local neighborhood density is calculated.
* First, the user can decide how the function "counts" other plants inside the buffer zone around the focal individual. Option 1) The function will tally the number of genets inside the buffer zone. Option 2) The function will calculate the proportion of the buffer zone that is occupied by other plants.
* Second, the user can decide whether the function will calculate intraspecific local neighborhood density (only consider other plants in the buffer zone of the same species as the focal individual) or interspecific local neighborhood density (consider all other plants in the buffer zone, regardless of species).
* Third, the user can decide whether the neighborhood density value (either a count or area) for each focal individual is a single value that sums the number or area of neighbors, or whether it is actually a list that provides the neighborhood density for each species present in the neighborhood.
```{r echo = FALSE, fig.width=6, fig.align = 'center', fig.cap = "Figure 3.1: This individual outlined in pink is a focal individual, and the pale pink shows a 10 cm buffer around it."}
datComp <- datTrackSpp[datTrackSpp$Year == 1922,]
ggplot(datComp) +
geom_sf(data =
sf::st_buffer(datComp[datComp$trackID == "HETCON_1922_6",], .1),
aes(), color = "#CC79A7", fill = "#CC79A7", alpha = .3, lty = 2) +
geom_sf(aes(fill = Species)) +
geom_sf(data = datComp[datComp$trackID == "HETCON_1922_6",],
aes(), color = "#CC79A7", fill = "#CC79A7", alpha = 0, lwd = 1.5) +
xlim(c(0,.5)) +
ylim(c(.2,.7)) +
xlab("quadrat horizontal edge (m)") +
ylab("quadrat vertical edge (m)") +
theme_classic() +
theme(legend.text = element_text(face = "italic"),
plot.margin = margin(1,0,1,0),
legend.title = element_blank()) +
scale_fill_discrete(type = c( "#009E73", "#F0E442", "#E69F00", "#56B4E9","#A6A6A6")) +
scale_color_discrete(type = c( "#009E73", "#F0E442","#E69F00", "#56B4E9", "#A6A6A6"))
```
```{r echo = FALSE, warnings = FALSE, include=FALSE }
suppressWarnings(st_intersects(st_buffer(datComp[datComp$trackID == "HETCON_1922_6",], .1), datComp))
countMethod <- datComp[c(3,5,63,67,70),]
areaMethodTemp <- suppressWarnings(suppressMessages(st_intersection(st_buffer(
datComp[datComp$trackID == "HETCON_1922_6",], .1), datComp)))
areaMethod <- areaMethodTemp[st_is(areaMethodTemp, c("POLYGON", "MULTIPOLYGON")),]
```
```{r echo = FALSE, warnings = FALSE, fig.width=6, fig.align = 'center', fig.cap = "Figure 3.2: The 10cm buffer around the focal individual overlaps with 5 other unique individuals of two species. These overlapping individuals are outlined in dark grey. Using the 'count' method in `getNeighbors()`, we would get an intraspecific competition value of 3, and an interspecific competition value of 5."}
ggplot(datComp) +
geom_sf(data =
suppressWarnings(sf::st_buffer(datComp[datComp$trackID == "HETCON_1922_6",], .1)),
aes(), color = "#CC79A7", fill = "#CC79A7", alpha = .3, lty = 2) +
geom_sf(aes(fill = Species)) +
geom_sf(data = datComp[datComp$trackID == "HETCON_1922_6",],
aes(), color = "#CC79A7", fill = "#CC79A7", alpha = 0, lwd = 1.5) +
geom_sf(data = countMethod, aes(), alpha = 0, lwd = 1.5) +
geom_sf_label(data = countMethod, label = c(1,2,3,4,5), nudge_x = .035,
label.padding = unit(.15, "lines"),
label.size = unit(.05,"mm")) +
xlim(c(0,.5)) +
ylim(c(.2,.7)) +
xlab("quadrat horizontal edge (m)") +
ylab("quadrat vertical edge (m)") +
theme_classic() +
theme(legend.text = element_text(face = "italic"),
plot.margin = margin(1,0,1,0),
legend.title = element_blank()) +
scale_fill_discrete(type = c( "#009E73", "#F0E442", "#E69F00", "#56B4E9","#A6A6A6")) +
scale_color_discrete(type = c( "#009E73", "#F0E442","#E69F00", "#56B4E9", "#A6A6A6"))
```
```{r echo = FALSE, warnings = FALSE, fig.width=6, fig.align = 'center', fig.cap = "Figure 3.3: The 10cm buffer around the focal individual overlaps with 5 other unique individuals of two species. The overlapping area is shaded in grey. Using the 'area' method in `getNeighbors()`, we would get an intraspecific competition metric of 0.0454, and an interspecific competition metric of 0.0462."}
areaMethodTemp <- suppressWarnings(st_intersection(st_buffer(
datComp[datComp$trackID == "HETCON_1922_6",], .1), datComp))
areaMethod <- areaMethodTemp[areaMethodTemp$trackID.1 != "HETCON_1922_6",]
ggplot(datComp) +
geom_sf(data =
suppressWarnings(sf::st_buffer(datComp[datComp$trackID == "HETCON_1922_6",], .1)),
aes(), color = "#CC79A7", fill = "#CC79A7", alpha = .3, lty = 2) +
geom_sf(aes(fill = Species)) +
geom_sf(data = datComp[datComp$trackID == "HETCON_1922_6",],
aes(), color = "#CC79A7", fill = "#CC79A7", alpha = 0, lwd = 1.5) +
geom_sf(data = areaMethod, aes(), fill = "#A6A6A6", color = "grey50", alpha = 0.7, lwd = 1.5) +
xlim(c(0,.5)) +
ylim(c(.2,.7)) +
xlab("quadrat horizontal edge (m)") +
ylab("quadrat vertical edge (m)") +
theme_classic() +
theme(legend.text = element_text(face = "italic"),
plot.margin = margin(1,0,1,0),
legend.title = element_blank()) +
scale_fill_discrete(type = c( "#009E73", "#F0E442", "#E69F00", "#56B4E9","#A6A6A6")) +
scale_color_discrete(type = c( "#009E73", "#F0E442","#E69F00", "#56B4E9", "#A6A6A6"))
## basal area of focal individual = 0.004981431
## area of buffer 0.06953875 - 0.004981431 = 0.06455732
## area of interspecific competitors 2.504885e-05 + 9.549092e-04 + 1.952117e-03 + 2.461883e-05 + 2.461883e-05 = 0.002981313
## area of intraspecific competitors 2.504885e-05 + 9.549092e-04 + 1.952117e-03 = 0.002932075
## proportion of interspecific competitors 0.002981313 / 0.06455732 = 0.04618087
## proportion of intraspecific competitors 0.002932075 / 0.06455732 = 0.04541816
```
Below are the arguments in the `getNeighbors()` function.
* **`dat`** An `sf` data frame in which each row represents data for a unique individual organism in a unique year. The `sf` geometry for each row must be either `MULTIPOLYGON` or `POLYGON` geometry. In addition to a "geometry" column, this data frame must have columns that contain data indicating the site, quadrat, site, and year of each observation. There also must be a column that contains a unique identifying value for each genet in each year. If `dat` is coming directly from `trackSpp()`, this column will be called "trackID".
* **`buff`** This is a single numeric value that indicates the desired width of the "buffer" around the focal individual in which the competitors are to be counted. This value must be in the same units as the spatial information in `dat`.
* **`method`** This is a character string that must equal either `"count"` or `"area"`. If `method = "count"`, then the number of individuals in the buffer area will be tallied. If `method = "area"`, then the proportion of the buffer area that is occupied by other individuals will be calculated.
* **`compType`** This is a character string that must be either `"allSpp"` or `"oneSpp"`. If `compType = "allSpp"`, then a metric of interspecific competition is calculated, meaning that every individual within the buffer around the focal individual is considered, no matter the species. If `compType = "oneSpp"`, then a metric of intraspecific competition is calculated, meaning that only individuals of the same species as the focal individual will be considered when calculating the competition metric. If no value is provided, it will default to "allSpp".
* **`output`** This is a character string that is set to either `"summed"` or `"bySpecies"`. The default is `"summed"`. This argument is only important to consider if you are using `compType = "allSpp"`. If `output = "summed"`, then only one count/area value is returned for each individual. This value is the total count or area of all neighbors within the focal species buffer zone, regardless of species. If `output = "bySpecies"`, there is a count or area value returned for each species present in the buffer zone. For example, you are using `getNeighbors()` with `method = "count"` and `compType = "allSpp"`. A focal individual in your dataset has seven other plants inside its buffer zone, three of species A, two of species B, and 2 of species C. If `output = "summed"`, the value in the "neighbors_count" column of the returned data frame will simply contain the value "7". If `output = "bySpecies"`, the "neighbors_count" column for this individual will actually contain a named list `{r} list("Species A "= 5, "Species B" = 3, "Species C" = 7)`. The default value of `output` is `"summed"`.
* **`trackID/species/quad/year/site/geometry`** These arguments only need to be included if the columns in `dat` that contain the data for trackID, species, site, quadrat, year and geometry of each observation are *different* from the names "trackID, "Species", "Site", "Quad", "Year, and "geometry". For example, if the column in your version of `dat` that contains the species identity of each observation is called "species_names", then the argument `species = "species_names"` must be included in your call to `getNeighbors()`.
<a id="getNeigh_out"></a>
#### *3.2* Function outputs
The output of `getNeighbors()` is an `sf` data frame that is identical to the input `dat`, but with either one or two additional columns. If `method = "area"`, there are two columns added called "nBuff_area" and "neighbors_area". The first contains the area of the buffer zone around each focal individual. The second contains the basal area of neighbors that overlap with a focal individual's buffer zone. If `method = "count"`, there is only one additional column added to the output, called "neighbors_count." This column contains a count of the neighbors that occur within a focal individual's buffer zone.
Here's an example of a `getNeighbors()` function call using the resulting data
from the example in section *2.2*, as well as the resulting data frame. Note
that `method = "area"`, so two columns are added to the returned data frame:
```{r}
datNeighbors <- plantTracker::getNeighbors(dat = datTrackSpp,
buff = .15,
method = "area",
compType = "allSpp")
```
```{r, echo = FALSE}
temp <- datNeighbors
rownames(temp) <- NULL
names(temp)[c(8,10,12)] <- c("basalArea", "survives_t+1", "size_t+1")
(temp)
#knitr::kable(head(temp[100:110,]), row.names = FALSE)
```
The example above uses the option `output = "summed"`, which is the default for the `getNeighbors()` function. With this option, the "neighbors_area" or "neighbors_count" column (depending on the `method` argument) contains just one value that sums the neighbor count or area across all neighbor species. However, if `output = "bySpecies"`, the "neighbors_count" or "neighbors_area" column contains a list with the counts or areas broken down by species. The `output` argument is described in more detail in section *3.1*. If you want to use the `getNeighbors()` function to determine how the effect of neighbors differs according to the species identity of those neighbors, setting `output = "bySpecies"` allows you to do this. However, it is likely that your subsequent analyses will need the by-species neighbor data in a matrix or data frame format, rather than in a list nested inside a data frame, which is how `getNeighbors()` returns the data.
Here is some code that turns the data returned by `getNeighbors(output = "bySpecies")` into a matrix where each row has data for one focal individual, and each column has data for one species. This format should be easier to work with for further analysis.
```{r}
# save the output of the getNeighbors() function
datNeighbors_bySpp <- plantTracker::getNeighbors(dat = datTrackSpp,
buff = .15, method = "area", compType = "allSpp", output = "bySpecies")
# determine all of the possible species that can occupy the buffer zone
compSpp <- unique(datTrackSpp$Species)
temp <- lapply(X = datNeighbors_bySpp$neighbors_area, FUN = function(x) {
tmp <- unlist(x)
tmp2 <- tmp[compSpp]
}
)
for (i in 1:length(temp)) {
# fix the column names
names(temp[[i]]) <- compSpp
# save the data in a matrix
if (i == 1) {
datOut <- temp[[i]]
} else {
datOut <- rbind(datOut, temp[[i]])
}
}
# make the rownames of the matrix correspond to the trackID of the focal individual
rownames(datOut) <- datNeighbors_bySpp$trackID
# show the first few rows of the datOut data frame:
datOut[1:5,]
```
___
<a id="next"></a>
### *4* Next steps
At this point, this dataset should be ready for you to use in any applications wish! There are a few additional functions that may help you in your analyses, and these are outlined in this section.
* Calculate recruitment by species-by-plot-by-year: the `getRecruits()` function
* Calculate population rate of increase (lambda) for each species in a plot: the `getLambda()` function
* Calculate the basal area of each species in each quadrat in each year: the `getBasalAreas()` function
* Make maps of quadrats over time: the `drawQuadMap()` function
Specific instructions for how to use each of these functions can be found in their documentation!