irg.R 3.87 KB
Newer Older
Alec L. Robitaille's avatar
Alec L. Robitaille committed
1 2
#' IRG
#'
3
#' Calculate the instantaneous rate of green-up.
Alec L. Robitaille's avatar
Alec L. Robitaille committed
4
#'
5
#' The DT argument expects a data.table of model estimated parameters for double logistic function of NDVI for each year and individual. Since it is the rate of green-up, model parameters required are only xmidS and scalS.
Alec L. Robitaille's avatar
Alec L. Robitaille committed
6
#'
7 8 9
#' The scaled argument is used to optionally rescale the IRG result to 0-1, for each year and individual.
#'
#' The id argument is used to split between sampling units. This may be a point id, polygon id, pixel id, etc. depending on your analysis. This should match the id provided to filtering functions.
10 11 12 13 14
#' The formula used is described in Bischoff et al. (2012):
#'
#' \deqn{IRG = (exp((t + xmidS) / scalS)) / (2 * scalS * (exp(1) ^ ((t + xmidS) / scalS)) + (scalS * (exp(1) ^ ((2 * t) / scalS))) + (scalS * exp(1) ^ ((2 * xmidS) / scalS)))}
#'
#' (See the "Getting started with irg vignette" for a better formatted formula.)
Alec L. Robitaille's avatar
Alec L. Robitaille committed
15 16 17 18 19 20 21
#'
#' @inheritParams model_ndvi
#' @inheritParams model_params
#' @param scaled boolean indicating if irg should be rescaled between 0-1 within id and year. If TRUE, provide id and year. Default is TRUE.
#'
#' @return
#'
22
#' Extended data.table 'irg' column of instantaneous rate of green-up calculated for each day of the year, for each individual and year.
Alec L. Robitaille's avatar
Alec L. Robitaille committed
23 24 25
#'
#' @export
#'
26 27
#' @family irg
#'
Alec L. Robitaille's avatar
Alec L. Robitaille committed
28
#' @examples
29 30 31 32 33 34 35 36 37 38 39
#' # Load data.table
#' library(data.table)
#'
#' # Read in example data
#' ndvi <- fread(system.file("extdata", "ndvi.csv", package = "irg"))
#'
#' # Filter and scale NDVI time series
#' filter_ndvi(ndvi)
#' scale_doy(ndvi)
#' scale_ndvi(ndvi)
#'
40 41 42
#' # Guess starting parameters
#' model_start(ndvi)
#'
43 44 45
#' # Double logistic model parameters given starting parameters for nls
#' mods <- model_params(
#'   ndvi,
46
#'   return = 'models',
Alec L. Robitaille's avatar
Alec L. Robitaille committed
47 48
#'   xmidS = 'xmidS_start',
#'   xmidA = 'xmidA_start',
49 50 51 52
#'   scalS = 0.05,
#'   scalA = 0.01
#' )
#'
53 54
#' # Fit double logistic curve to NDVI time series
#' fit <- model_ndvi(mods, observed = FALSE)
Alec L. Robitaille's avatar
Alec L. Robitaille committed
55 56
#'
#' # Calculate IRG for each day of the year
57
#' calc_irg(fit)
Alec L. Robitaille's avatar
fst irg  
Alec L. Robitaille committed
58 59 60
calc_irg <-
	function(DT,
					 id = 'id',
61 62
					 year = 'yr',
					 scaled = TRUE) {
Alec L. Robitaille's avatar
fst irg  
Alec L. Robitaille committed
63 64
		# NSE error
		xmidS <- scalS <- irg <- NULL
Alec L. Robitaille's avatar
Alec L. Robitaille committed
65

66
		check_truelength(DT)
Alec L. Robitaille's avatar
fst irg  
Alec L. Robitaille committed
67 68
		check_col(DT, 'xmidS')
		check_col(DT, 'scalS')
69
		check_col(DT, 't')
70

Alec L. Robitaille's avatar
fst irg  
Alec L. Robitaille committed
71 72 73
		if (any(unlist(DT[, lapply(.SD, function(x)
			any(is.na(x)))]))) {
			warning('NAs found in DT, IRG will be set to NA.')
Alec L. Robitaille's avatar
Alec L. Robitaille committed
74 75
		}

Alec L. Robitaille's avatar
fst irg  
Alec L. Robitaille committed
76 77 78 79 80 81 82
		DT[, irg :=
			 	(exp((t + xmidS) / scalS)) /
			 	(2 * scalS * (exp(1) ^ ((t + xmidS) / scalS)) +
			 	 	(scalS * (exp(1) ^ ((
			 	 		2 * t
			 	 	) / scalS))) +
			 	 	(scalS * exp(1) ^ ((2 * xmidS) / scalS)))]
83

Alec L. Robitaille's avatar
fst irg  
Alec L. Robitaille committed
84 85 86
		if (scaled) {
			check_col(DT, id, 'id')
			check_col(DT, year, 'year')
Alec L. Robitaille's avatar
Alec L. Robitaille committed
87

Alec L. Robitaille's avatar
Alec L. Robitaille committed
88 89 90
			DT[, irg :=
				 	(irg - min(irg, na.rm = TRUE)) /
				 	(max(irg, na.rm = TRUE) - min(irg, na.rm = TRUE)),
Alec L. Robitaille's avatar
fst irg  
Alec L. Robitaille committed
91 92
				 by = c(id, year)]
		}
Alec L. Robitaille's avatar
Alec L. Robitaille committed
93

Alec L. Robitaille's avatar
fst irg  
Alec L. Robitaille committed
94
		return(DT)
Alec L. Robitaille's avatar
Alec L. Robitaille committed
95 96
	}

Alec L. Robitaille's avatar
fst irg  
Alec L. Robitaille committed
97 98 99 100


#' IRG
#'
Alec L. Robitaille's avatar
Alec L. Robitaille committed
101
#' Wrapper function for one step IRG calculation. Only defaults.
Alec L. Robitaille's avatar
fst irg  
Alec L. Robitaille committed
102
#'
Alec L. Robitaille's avatar
Alec L. Robitaille committed
103
#' data.table must have columns:
104
#'
Alec L. Robitaille's avatar
Alec L. Robitaille committed
105 106 107 108 109 110 111
#' \itemize{
#'   \item 'id' - individual identifier
#'   \item 'yr' - year of observation
#'   \item 'NDVI' - NDVI value
#'   \item 'DayOfYear' - day of year/julian day of observation
#'   \item 'SummaryQA' - summary quality value for each sample (provided by MODIS)
#' }
112 113 114
#'
#' @inheritParams filter_qa
#'
Alec L. Robitaille's avatar
fst irg  
Alec L. Robitaille committed
115 116 117 118 119 120 121 122 123 124 125 126 127 128 129
#' @return
#'
#' Extended data.table 'irg' column of instantaneous rate of green-up calculated for each day of the year, for each individual and year.
#'
#' @export
#'
#' @family irg
#'
#' @examples
#' # Load data.table
#' library(data.table)
#'
#' # Read in example data
#' ndvi <- fread(system.file("extdata", "ndvi.csv", package = "irg"))
#'
Alec L. Robitaille's avatar
Alec L. Robitaille committed
130
#' # Calculate IRG for each day of the year and individual
Alec L. Robitaille's avatar
fst irg  
Alec L. Robitaille committed
131 132
#' out <- irg(ndvi)
irg <- function(DT) {
Alec L. Robitaille's avatar
Alec L. Robitaille committed
133 134
	check_truelength(DT)

135 136 137 138
	filter_ndvi(DT)
	scale_doy(DT)
	scale_ndvi(DT)
	model_start(DT)
Alec L. Robitaille's avatar
Alec L. Robitaille committed
139
	model_params(DT, returns = 'columns',
Alec L. Robitaille's avatar
Alec L. Robitaille committed
140
							 xmidS = 'xmidS_start', xmidA = 'xmidA_start',
Alec L. Robitaille's avatar
Alec L. Robitaille committed
141 142
							 scalS = 0.05, scalA = 0.1)
	return(calc_irg(model_ndvi(DT, observed = FALSE)))
Alec L. Robitaille's avatar
Alec L. Robitaille committed
143
}