Prototyping with tools    Next: Generating Normal(0,1) Random Variables Up: Octave Programs Previous: Octave Programs   Index

Click for printer friendely version of this HowTo

# Univariate General Linear Models

```function beta=general_linear(x,y,c)

# y is a column (nx1) vector of observations
# x is the design matrix (nxm)
# c is the contrast matrix
# NOTE: it is assumed that c is set up so that theta = 0

n = max(size(x)); # number of observations
p = min(size(x)); # number of parameters
num_tests = min(size(c)); # number of tests

beta = inv(x'*x)*x'*y # estimate the parameters

# test model
numerator = beta'*c'*inv((c*inv(x'*x)*c'))*c*beta / num_tests;
denominator = (y'*y - y'*x*beta)/(n-p); # this is also known as MSE
f_test = numerator / denominator
p_value = 1 - f_cdf(f_test, num_tests, n-p) # calculate p-value

if (num_tests == 1)
t_test = sqrt(f_test)
end

# generate graphical output
MIN_F_DIST = 0;
MAX_F_DIST = 20;
NUM_STEPS = 100;

graph_x = linspace(MIN_F_DIST, MAX_F_DIST, NUM_STEPS);
graph_y = f_pdf(graph_x, num_tests, n-p);

description = sprintf("g;F(%d, %d)-Dist;", num_tests, n-p);
plot(graph_x, graph_y, description, f_test, 0, "r*;Your Data;")

endfunction
```

Frank Starmer 2004-05-19
>