forked from egarpor/goffda
-
Notifications
You must be signed in to change notification settings - Fork 0
/
Copy pathREADME.Rmd
136 lines (104 loc) · 4.96 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
---
output:
md_document:
variant: gfm
---
```{r, echo = FALSE}
knitr::opts_chunk$set(
collapse = TRUE, comment = "#>", fig.path = "README/README-",
message = FALSE, warning = FALSE, fig.asp = 1, fig.align = 'center'
)
```
goffda
======
<!-- badges: start -->
[](https://www.gnu.org/licenses/gpl-3.0)
[](https://travis-ci.org/egarpor/goffda)
[](https://cran.r-project.org/package=goffda)
<!-- badges: end -->
## Overview
Software companion for the paper "*A goodness-of-fit test for the functional linear model with functional response*" (García-Portugués, Álvarez-Liébana, Álvarez-Pérez and González-Manteiga, 2019). It implements the proposed estimators and goodness-of-fit tests for the functional linear model with scalar response. It also allows to replicate the data application presented.
## Installation
```{r, install, eval = FALSE}
# Install the package
library(devtools)
install_github("egarpor/goffda")
# Load package
library(goffda)
```
```{r, load, echo = FALSE}
# Load package
library(goffda)
```
## Usage
The following are some simple examples of the usage of the main function of the package, `flm_test`. More examples are available in `?flm_test`.
```{r, test}
# Generate data under H0
n <- 100
set.seed(987654321)
X_fdata <- r_ou(n = n, t = seq(0, 1, l = 101), sigma = 2)
epsilon <- r_ou(n = n, t = seq(0, 1, l = 101), sigma = 0.5)
Y_fdata <- epsilon
# Test the FLMFR
flm_test(X = X_fdata, Y = Y_fdata, verbose = FALSE)
# Simple hypothesis
flm_test(X = X_fdata, Y = Y_fdata, beta0 = 0, verbose = FALSE)
# Generate data under H1
n <- 100
set.seed(987654321)
sample_frm_fr <- r_frm_fr(n = n, scenario = 3, s = seq(0, 1, l = 101),
t = seq(0, 1, l = 101), nonlinear = "quadratic")
X_fdata <- sample_frm_fr[["X_fdata"]]
Y_fdata <- sample_frm_fr[["Y_fdata"]]
# Test the FLMFR
flm_test(X = X_fdata, Y = Y_fdata, verbose = FALSE)
```
## AEMET temperatures data
An illustration of the test application in a real dataset.
```{r, aemet}
## Data preprocessing
# Load raw data
data("aemet_temp")
# Partition the dataset
with(aemet_temp, {
ind_pred <- which((1974 <= df$year) & (df$year <= 1993))
ind_resp <- which((1994 <= df$year) & (df$year <= 2013))
aemet_temp_pred <<- list("df" = df[ind_pred, ], "temp" = temp[ind_pred])
aemet_temp_resp <<- list("df" = df[ind_resp, ], "temp" = temp[ind_resp])
})
# Average the temperature on each period
mean_aemet <- function(x) {
m <- tapply(X = 1:nrow(x$temp$data), INDEX = x$df$ind,
FUN = function(i) colMeans(x$temp$data[i, , drop = FALSE],
na.rm = TRUE))
x$temp$data <- do.call(rbind, m)
return(x$temp)
}
# Build predictor and response fdatas
aemet_temp_pred <- mean_aemet(aemet_temp_pred)
aemet_temp_resp <- mean_aemet(aemet_temp_resp)
# Average yearly temperatures
avg_year_pred <- rowMeans(aemet_temp_pred$data)
avg_year_resp <- rowMeans(aemet_temp_resp$data)
# Average temperatures on both periods
plot(func_mean(aemet_temp_pred), ylim = c(5, 30),
main = "Average AEMET temperature")
plot(func_mean(aemet_temp_resp), col = 2, add = TRUE)
legend("topleft", legend = c("1974-1993", "1994-2013"), col = 1:2, lwd = 2)
# Test composite and simple hypothesis
(gof <- flm_test(X = aemet_temp_pred, Y = aemet_temp_resp, B = 1e4,
verbose = FALSE, plot_dens = FALSE, plot_proc = FALSE))
beta0 <- diag(x = rep(1, length(aemet_temp_pred$argvals)))
flm_test(X = aemet_temp_pred, Y = aemet_temp_resp, verbose = FALSE,
beta0 = beta0, B = 1e4, plot_dens = FALSE, plot_proc = FALSE)
# Visualize estimation
lev <- seq(-0.022, 0.022, by = 0.004)
col <- colorRampPalette(colors = c("blue", "white", "red"))
filled.contour(x = aemet_temp_pred$argvals, y = aemet_temp_resp$argvals,
z = gof$fit_flm$Beta_hat, levels = lev, color.palette = col)
```
## References
García-Portugués, E., Álvarez-Liébana, J., Álvarez-Pérez, G. and González-Manteiga, W. (2019). A goodness-of-fit test for the functional linear model with functional response. *arXiv:1909.07686*. <https://arxiv.org/abs/1909.07686>
Cuesta-Albertos, J. A., García-Portugués, E., Febrero-Bande, M. and González-Manteiga, W. (2019). Goodness-of-fit tests for the functional linear model based on randomly projected empirical processes. *Annals of Statistics*, 47(1):439-467. <https://doi.org/10.1214/18-AOS1693>
Febrero-Bande, M. and Oviedo de la Fuente, M. (2012). Statistical Computing in Functional Data Analysis: The R Package fda.usc. *Journal of Statistical Software*, 51(4):1-28. <http://www.jstatsoft.org/v51/i04/>
García-Portugués, E., González-Manteiga, W. and Febrero-Bande, M. (2014). A goodness-of-fit test for the functional linear model with scalar response. *Journal of Computational and Graphical Statistics*, 23(3):761-778. <http://doi.org/10.1080/10618600.2013.812519>