# R/calcs.r In growthPheno: Functional Analysis of Phenotypic Growth Data to Smooth and Extract Traits

```#Various functions for doing calculations.

#Function to calculate the cumulative sum, ignoring the first element if exclude.1st is TRUE
"cumulate" <- function(x, exclude.1st=FALSE, na.rm = FALSE, ...)
{
if (na.rm)
x <- x[!(is.na(x))]
sum <- x
if (exclude.1st)
sum[-1] <- cumsum(x[-1])
else
sum <- cumsum(x)
return(sum)
}

#Functions to calculate growth rates between successive imagings
"AGRdiff" <- function(x, time.diffs, lag=1)
{
x.diffs <- calcLagged(x, operation = "-", lag = lag)
x.diffs <- x.diffs / time.diffs
}
"PGR" <- function(x, time.diffs, lag=1)
{
x.rates <- calcLagged(x, operation = "/", lag = lag)
x.rates <- x.rates ^ (1/time.diffs)
}

"RGRdiff" <- function(x, time.diffs, lag=1)
{
x.rates <- log(PGR(x, time.diffs, lag = lag))
}

"WUI" <- function(response, water)
{
response.WUI <- ifelse(water != 0,
response / water,
NA)
return(response.WUI)
}

#Functions to calculate a single-valued function, including the observation has the value of the function

"which.funct.value" <- function(x, FUNCT = NULL, na.rm = TRUE, ...)
{
funct <- get(FUNCT)
funct <- match.fun(funct)
#  k <- which(x==funct(x, ...))
if (all(is.na(x)))
k <- NA
else
{
adiff <- abs(x - funct(x, na.rm = na.rm, ...))
if (length(k) == 0)
k <- NA
else if (length(k) > 1)
k <- k[1]
}
return(k)
}
#Functions to calculate statistics for a profile

#Function to calculate the sample range
#"sample.range" <- function(x, na.rm = FALSE){z <- diff(range(x, na.rm = na.rm))}

#Function to calculate the root mean square
#"rms" <- function(x, na.rm = FALSE){z <- sqrt(sum(x*x, na.rm = na.rm)/length(x))}

#Function to calculate the sum over an interval between start.value and end.value
"interval.sum" <- function(x, start.value=1, end.value=1, na.rm = FALSE)
# Get the sum of x for which subset.var has values between start.values and end.values
{ z <- sum(x[start.value:end.value], na.rm = na.rm)
}

#Function to calculate the sum over an interval between start.value and end.value
"interval.wtmean" <- function(x, w, start.value=1, end.value=1, na.rm = FALSE)
# Get the sum of x for which subset.var has values between start.values and end.values
{ z <- sum(w[start.value:end.value] * x[start.value:end.value], na.rm = na.rm) /
sum(w[start.value:end.value])
}

#Function to get the values for specified times and, optionally, a column with the times factor
"getTimesSubset" <- function(data, responses,
individuals = "Snapshot.ID.Tag", times = "DAP",
which.times, suffix = NULL, sep.suffix.times = ".",
include.times = FALSE, include.individuals = FALSE)
{
n <- dim(data)[1]
subset <- rep(FALSE, n)
for (day in which.times)
{
this.times <- data[times] == day
subset <- subset | this.times
}
data.sub <- data[subset, ]
data.sub <- data.sub[do.call(order, data.sub), ]
if (!is.null(suffix))
new.responses <- unlist(lapply(responses,
function(name, suffix){paste(name, suffix, sep=sep.suffix.times)},
suffix=suffix))
else
new.responses <- responses
if (include.individuals)
{
responses <- c(responses,individuals)
new.responses <- c(new.responses,individuals)
}
if (include.times)
{
responses <- c(responses,times)
if (!is.null(suffix))
new.responses <- c(new.responses, paste(times, suffix, sep = sep.suffix.times))
else
new.responses <- c(new.responses,times)
}
data.sub <- data.sub[responses]
names(data.sub) <- new.responses
return(data.sub)
}

#Functions to do calculations between successive dates
# - does not assume same number time points for all individuals
#"Replace"  <- function(x, y) {z <- y}
"calcLagged" <- function(x, operation = NULL, lag=1)
#This function replaces the observations with values calculated
# (i) for positive lag, itself and the value lag observations before it,
# (ii) for negative lag, itself and the value lag observations after it.
#operation specifies calculation to be made on the pair of  values
#It returns as many values as are in data, the 1st lag values being NA
{
n <- length(x)
nl <- n-abs(lag)
if (is.null(operation))
{
if (lag > 0)
x[(lag+1):n] <- x[1:nl]
else
{
if (lag < 0)
x[1:nl] <- x[(abs(lag)+1):n]
}
}
else
{
FUN <- get(operation)
FUN <- match.fun(FUN)
if (lag > 0)
x[(lag+1):n] <- FUN(x[(lag+1):n], x[1:nl])
else
{
if (lag < 0)
x[1:nl] <- FUN(x[1:nl], x[(abs(lag)+1):n])
else
x[1:n] <- FUN(x[1:n], x[1:n])
}
}
if (lag > 0)
x[1:lag] <- NA
else
x[(nl+1):n] <- NA
return(x)
}

#Function to test if any values in a set of values are anomalous
# in being outside specified limits
"anom" <- function(x, lower=NULL, upper=NULL, na.rm = TRUE)
{ if (is.null(lower))
{ if (is.null(upper))
stop("Must supply at least a lower or an  upper limit below or above which values are anomalous")
else
anom <- any(x > upper, na.rm=na.rm)
} else
{ if (is.null(upper))
anom <- any(x < lower, na.rm=na.rm)
else
anom <- any(x < lower, na.rm=na.rm) || any(x > upper, na.rm=na.rm)
}
return(anom)
}

"twoLevelOpcreate" <- function(data, responses, treatment.factor="Treatment.1",
suffices.treatment=c("Cont","Salt"), control=1,
columns.suffixed = NULL,
operations = "/", suffices.results="OST",
columns.retained = c("Snapshot.ID.Tag",
"Smarthouse","Lane",
"Zone","cZone","SHZone",
"ZLane","ZMainunit","cMainPosn",
"Genotype.ID"),
by = c("Smarthouse","Zone","ZMainunit"))
{
#Check that response and treatment.factor are in data
checkNamesInData(c(responses, treatment.factor), data)

trt.lev <- levels(data[[treatment.factor]])
ntrt <- length(trt.lev)
if (ntrt != 2)
stop("The treatment factor must have only two levels")
if (!(control %in% c(1,2)))
stop("contol must be either 1 or 2")
if (!all(by %in% columns.retained))
stop("The names in by must be included in columns.retained")
if (treatment.factor %in% columns.retained)
stop("The treatment.factor should not be included in columns.retained")
treated <- c(1,2)[c(1,2)!=control]
if (!is.null(columns.suffixed))
if (any(!(unique(columns.suffixed) %in% c(responses,columns.retained))))
stop("The columns.suffixed must be included in columns.retained")

#identify the functions specified in operations and check the number of them
nresp <-  length(responses)
fun <- vector(mode="list", length=nresp)
fun <- lapply(operations,
function(operations)
fun <- match.fun(get(operations))
)
if (length(operations)==1 && nresp != 1)
fun[2:nresp] <- fun[1]
if (length(fun) != nresp)
stop("The length of operations is neither 1 nor the length of responses")
if (length(suffices.results) == 1 && nresp != 1)
suffices.results <- rep(suffices.results, nresp)
if (length(suffices.results) != nresp)
stop("The length of suffices.results is neither 1 nor the length of responses")

#Split data according to the two treatments
uresponses <- unique(responses)
ucols.suffixed <- unique(columns.suffixed)
data <- data[c(treatment.factor,columns.retained, uresponses)]
data <- split(data, data[[treatment.factor]])
data[[1]] <- data[[1]][-match(treatment.factor,names(data[[1]]))]
data[[2]] <- data[[2]][-match(treatment.factor,names(data[[2]]))]

#Create data frame with the retained columns and columns for the GR values for each treatment
data.out <- data[[treated]]
names(data.out)[match(c(uresponses,ucols.suffixed), names(data.out))] <-
paste(c(uresponses,ucols.suffixed),suffices.treatment[treated],sep=".")
data.out <- merge(data.out, data[[control]][c(by,uresponses,ucols.suffixed)], by=by,sort=FALSE)
names(data.out)[match(c(uresponses,ucols.suffixed), names(data.out))] <-
paste(c(uresponses,ucols.suffixed),suffices.treatment[control],sep=".")

#Add the results of the operations
data.out[paste(responses,suffices.results,sep=".")] <-
lapply(1:nresp,
function(k, responses, fun, data)
{ out <- fun[[k]](data[[paste(responses[[k]],suffices.treatment[treated],sep=".")]],
data[[paste(responses[[k]],suffices.treatment[control],sep=".")]])
return(out)
}, responses=responses, fun=fun, data=data.out)
return(data.out)
}
```

## Try the growthPheno package in your browser

Any scripts or data that you put into this service are public.

growthPheno documentation built on May 29, 2024, 6:03 a.m.