ph_test_statistic {gmtFD}

Pointwise Hotelling's T^2-test statistic

Description

The function ph_test_statistic() calculates the pointwise Hotelling's T^2-test statistic.

Usage

ph_test_statistic(x, h)

Arguments

x

a list of k elements corresponding to groups. Each element representing a group is a list of p elements corresponding to functional variables, and each such element (representing a functional variable) is a matrix of size n_i\times ntp of descrete observations in design time points. ntp denotes a number of design time points.

h

contrast matrix. For contrast matrices based on Dunnett’s and Tukey’s contrasts, it can be created by the contr_mat() function from the package GFDmcv (see examples).

Details

For details, see the documentation of the gmtFD() function or the papers Munko et al. (2023, 2024).

Value

A vector of values of the pointwise Hotelling's T^2-test statistic.

References

Dunnett C. (1955) A multiple comparison procedure for comparing several treatments with a control. Journal of the American Statistical Association 50, 1096-1121.

Munko M., Ditzhaus M., Pauly M., Smaga L., Zhang J.T. (2023) General multiple tests for functional data. Preprint https://arxiv.org/abs/2306.15259

Munko M., Ditzhaus M., Pauly M., Smaga L. (2024) Multiple comparison procedures for simultaneous inference in functional MANOVA. Preprint https://arxiv.org/abs/2406.01242

Tukey J.W. (1953) The problem of multiple comparisons. Princeton University.

Examples

# Some of the examples may run some time.

# Canadian weather data set
# There are three samples of mean temperature and precipitations for
# fifteen weather stations in Western Canada, another fifteen in Eastern Canada, 
# and the remaining five in Northern Canada.

# one functional variable - temperature
library(fda)
data_set_t <- t(CanadianWeather$dailyAv[,, "Temperature.C"])
# number of samples
k <- 3
# number of variables
p <- 1
# preparing data set
gr_label <- rep(c(1, 2, 3), c(15, 15, 5))
data_set <- list(list(data_set_t[gr_label == 1, ]),
                 list(data_set_t[gr_label == 2, ]),
                 list(data_set_t[gr_label == 3, ]))

# Tukey's contrast matrix
h_tukey <- GFDmcv::contr_mat(k, type = "Tukey")
h_tukey_m <- kronecker(h_tukey, diag(p))
# plots for pointwise Hotelling's T^2-test statistics
oldpar <- par(mfrow = c(2, 2), mar = c(4, 2, 2, 0.1))
plot(ph_test_statistic(data_set, h_tukey_m), type = "l",
     ylim = c(0, max(ph_test_statistic(data_set, h_tukey_m))),
     main = "Global hypothesis", xlab = "Day")
plot(ph_test_statistic(data_set, matrix(h_tukey_m[1, ], 1)), type = "l",
     ylim = c(0, max(ph_test_statistic(data_set, h_tukey_m))),
     main = "Contrast 1", xlab = "Day")
plot(ph_test_statistic(data_set, matrix(h_tukey_m[2, ], 1)), type = "l",
     ylim = c(0, max(ph_test_statistic(data_set, h_tukey_m))),
     main = "Contrast 2", xlab = "Day")
plot(ph_test_statistic(data_set, matrix(h_tukey_m[3, ], 1)), type = "l",
     ylim = c(0, max(ph_test_statistic(data_set, h_tukey_m))),
     main = "Contrast 3", xlab = "Day")
par(oldpar)

# Dunnett's contrast matrix
h_dunnett <- GFDmcv::contr_mat(k, type = "Dunnett")
h_dunnett_m <- kronecker(h_dunnett, diag(p))
# plots for pointwise Hotelling's T^2-test statistics
oldpar <- par(mfrow = c(3, 1), mar = c(4, 2, 2, 0.1))
plot(ph_test_statistic(data_set, h_dunnett_m), type = "l",
     ylim = c(0, max(ph_test_statistic(data_set, h_dunnett_m))),
     main = "Global hypothesis", xlab = "Day")
plot(ph_test_statistic(data_set, matrix(h_dunnett_m[1, ], 1)), type = "l",
     ylim = c(0, max(ph_test_statistic(data_set, h_dunnett_m))),
     main = "Contrast 1", xlab = "Day")
plot(ph_test_statistic(data_set, matrix(h_dunnett_m[2, ], 1)), type = "l",
     ylim = c(0, max(ph_test_statistic(data_set, h_dunnett_m))),
     main = "Contrast 2", xlab = "Day")
par(oldpar)

# two functional variables - temperature and precipitation
library(fda)
data_set_t <- t(CanadianWeather$dailyAv[,, "Temperature.C"])
data_set_p <- t(CanadianWeather$dailyAv[,, "Precipitation.mm"])
# number of samples
k <- 3
# number of variables
p <- 2
# preparing data set
gr_label <- rep(c(1, 2, 3), c(15, 15, 5))
data_set <- list(list(data_set_t[gr_label == 1, ], data_set_p[gr_label == 1, ]),
                 list(data_set_t[gr_label == 2, ], data_set_p[gr_label == 2, ]),
                 list(data_set_t[gr_label == 3, ], data_set_p[gr_label == 3, ]))

# Tukey's contrast matrix
h_tukey <- GFDmcv::contr_mat(k, type = "Tukey")
h_tukey_m <- kronecker(h_tukey, diag(p))
# plots for pointwise Hotelling's T^2-test statistics
oldpar <- par(mfrow = c(2, 2), mar = c(4, 2, 2, 0.1))
plot(ph_test_statistic(data_set, h_tukey_m), type = "l",
     ylim = c(0, max(ph_test_statistic(data_set, h_tukey_m))),
     main = "Global hypothesis", xlab = "Day")
plot(ph_test_statistic(data_set, h_tukey_m[1:2, ]), type = "l",
     ylim = c(0, max(ph_test_statistic(data_set, h_tukey_m))),
     main = "Contrast 1", xlab = "Day")
plot(ph_test_statistic(data_set, h_tukey_m[3:4, ]), type = "l",
     ylim = c(0, max(ph_test_statistic(data_set, h_tukey_m))),
     main = "Contrast 2", xlab = "Day")
plot(ph_test_statistic(data_set, h_tukey_m[5:6, ]), type = "l",
     ylim = c(0, max(ph_test_statistic(data_set, h_tukey_m))),
     main = "Contrast 3", xlab = "Day")
par(oldpar)

# Dunnett's contrast matrix
h_dunnett <- GFDmcv::contr_mat(k, type = "Dunnett")
h_dunnett_m <- kronecker(h_dunnett, diag(p))
# plots for pointwise Hotelling's T^2-test statistics
oldpar <- par(mfrow = c(3, 1), mar = c(4, 2, 2, 0.1))
plot(ph_test_statistic(data_set, h_dunnett_m), type = "l",
     ylim = c(0, max(ph_test_statistic(data_set, h_dunnett_m))),
     main = "Global hypothesis", xlab = "Day")
plot(ph_test_statistic(data_set, matrix(h_dunnett_m[1, ], 1)), type = "l",
     ylim = c(0, max(ph_test_statistic(data_set, h_dunnett_m))),
     main = "Contrast 1", xlab = "Day")
plot(ph_test_statistic(data_set, matrix(h_dunnett_m[2, ], 1)), type = "l",
     ylim = c(0, max(ph_test_statistic(data_set, h_dunnett_m))),
     main = "Contrast 2", xlab = "Day")
par(oldpar)