Nothing
###################################################################
#
# These functions are part of WACSgen V1.0
# Copyright © 2013,2014,2015, D. Allard, BioSP,
# and Ronan Trépos MIA-T, INRA
#
# This program is free software; you can redistribute it and/or
# modify it under the terms of the GNU General Public License
# as published by the Free Software Foundation; either version 2
# of the License, or (at your option) any later version.
#
# This program is distributed in the hope that it will be useful,
# but WITHOUT ANY WARRANTY; without even the implied warranty of
# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
# GNU General Public License for more details. http://www.gnu.org
#
###################################################################
wacs.estimcycle = function(y,spar=0.7,trend.norm="L2",plot.it=TRUE,DIR="./"){
#####################################################
#
# wacs.estimcycle : removes the seasonal component of the values.
# returns the residuals, the central tendency
# and a measure of deviation
#
# ARGUMENT
# y : vector of values; length must be a multiple of 365
# spar : smoothness paramater for seasonal trend
# trend.norm: norm to be used
# "L2" = quadratic norm, leading to mean and st
# "L1" = Absolute value norm, leading to median and absolute deviation
# plot.it : whether plots are produced or not
# DIR : directory in which plots are printed
#
# VALUE
# A list containing three elements:
# 1/ an array containing the residuals (Nrow = Nd, number of days;
# Ncol = Nv, number of variables to be centered)
# 2/ 365 x Nv array of central tendencies.
# 3/ 365 x Nv array of central deviations.
#
# CALLS
# extract.annual.trend
#####################################################
Central = NULL
Deviation = NULL
Resid = NULL
Nv = dim(y)[2]
for (v in 1:Nv)
{
M = matrix(y[,v],ncol=365,byrow=TRUE)
Dv = extract.annual.trend(y[,v],spar=spar,trend.norm=trend.norm)
Central = cbind(Central,Dv[,1])
Deviation = cbind(Deviation,Dv[,2])
Resid = cbind(Resid,(y[,v] - Dv[,1])/Dv[,2])
colnames(Resid)[v] = paste("V",v,"d",sep="")
if (plot.it==TRUE){
grDevices::png(paste(DIR,"Annual_",colnames(y)[v],".png",sep=""))
graphics::boxplot(t(M)~c(1:365),outline=FALSE,xlab="Julian day",
ylab=paste("Variable V",v,sep=""),na.action=drop)
graphics::lines(Dv[,1],col="red",lwd=3)
graphics::lines(Dv[,1]-2*Dv[,2],col="red",lwd=1)
graphics::lines(Dv[,1]+2*Dv[,2],col="red",lwd=1)
grDevices::dev.off()
}
}
return(list(Resid,Central,Deviation))
}
extract.annual.trend = function(y,spar=0.7,trend.norm="L2"){
#####################################################
#
# extract.annual.trend :extract the central tendency and a measure of deviation
#
# ARGUMENTS
# y : vector of values; must be of length 365
# spar: smoothness paramater for seasonal trend
# trend.norm: norm to be used
# "L2" = quadratic norm, leading to mean and st
# "L1" = Absolute value norm, leading to median and absolute deviation
#
# VALUE
# a 365 x 2 array. First column: vector of central values;
# Second column: vector of deviation
#
#####################################################
M=matrix(y,ncol=365,byrow=TRUE)
if (trend.norm=="L1"){
y.center = apply(M,2,stats::median,na.rm=TRUE)
dev = function(v) z=mean(abs(v-stats::median(v,na.rm=TRUE)),na.rm=TRUE)
}else{}
if (trend.norm=="L2"){
y.center = apply(M,2,mean,na.rm=TRUE)
dev = function(v) z = sqrt(mean( (v-mean(v,na.rm=TRUE))^2, na.rm=TRUE))
}else{}
y.dev = apply(M,2,dev)
# smoothing
y.center.smooth = stats::smooth.spline(rep(y.center,3),spar=spar)[[2]][366:730]
y.dev.smooth = stats::smooth.spline(rep(y.dev,3),spar=spar)[[2]][366:730]
return(cbind(y.center.smooth,y.dev.smooth))
}
Any scripts or data that you put into this service are public.
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.