We can perform this by simply calculating the Hotelling T2 value for each observation of the data we generate and then comparing it to the threshold obtained from a chi-sqaure distribution. Then we can count the number of observations which exceed the threshold and determine the probability of shift detection using the result:
Here is the workflow with code
- Create mean vector and the identity correlation matrix for the variables to create data without shift. Then add shift to the first variable.
mu = zeros(1, v); % Mean vector
Sigma = eye(v); % Identity correlation matrix
data = mvnrnd(mu, Sigma, n); % Generate data without shift
data(:, 1) = data(:, 1) + shiftMagnitude; % Add shift to x1 variable
%Note:
% v is the number of variables
% n is the number of observations we have for each variable.
%To obtain clear results keep n > 1000
- Calculate the Hostelling T2 statistic for each observation.
T2 = zeros(n, 1);
for i = 1:n
x = data(i, :);
T2(i) = (x - mu) * inv(Sigma) * (x - mu)';
end
- Calculate the threshold from the chi-square distribution (keep the false probability = 0.002).
falseAlarm = 0.002; % False alarm probability
threshold = chi2inv(1 - falseAlarm, p);
- Calculate the number observations with T2 value greater than threshold.
detections = sum(T2 > threshold);
- Calculate the probability of shift detection.
probabilityOfDetection = detections/n;
