linstats 2006b

Michael Boedigheimer (view profile)

27 Dec 2006 (Updated )

linear multivariate statistics

tests2eqn( glm, tests, reference, short_form )
```function source = tests2eqn( glm, tests, reference, short_form )
% TESTS2EQN  produces names for anova tests (helper for anova)
%
% function source = tests2eqn( glm, tests, reference )
%
% In the function 'anova' tests are encoded in matrix form.
% For readability (and debugging) this function is provided to convert
% the matrices into a human readable form.
%
% An anova test compares two models. One model, the full model, is reduced
% compared to the other model, reduced model, in that it contains a subset
% of the terms in the reference model. The terms that are left out are the
% terms that are being tests.
% See gettests
%
% input
% each row of test is a matrix of terms that will dropped from the model. A
% value of 1 indicates that the corresponding term will be dropped, and
% this will become the reduced model. It will be compared to a reference
% model (full model), which is specified in another row in the tests
% matrix. Each reduced model is compared to a full model, which is another
% row in the tests matrix. The row is given in the vector 'reference'. Thus
% the reference vector specifies the number of tests to be done. If there
% are n tests  then the first n rows of the tests matrix specify the
% reduced models and the ith reduced model is compared to the row in
% reference(i). Thus, the two models in the ith anova test are tests(i,:)
% and tests(reference(i),:)
%
%   for example say there are three terms in the model, x1,x2,x3
%   x1|x2,x3     - test the effect (reduction in sse) by adding the term
%                  x1 to a model already containing x2 and x3
%   x1,x2|x3     - test the effect of adding x1 and x2 to a model already
%                  containing x3
%   x1|Intercept - test the effect of adding x1 to the null model
%   x1(x2)|x3    - test the effect of adding x1 to a model containing x3
%                  leaving out the effect of x2
%
% Example
% glm = encode( X, 3, 2, Origin, Cylinders );  % Cylinders is explanatory variable
% [tests, reference] = gettests( glm, 1 );
% tests2eqn(glm, tests, reference )
% [tests, reference] = gettests( glm, 3 );
% tests2eqn(glm, tests, reference )
% tests2eqn(glm, tests, reference, 1 )

% \$Id: tests2eqn.m,v 1.6 2006/12/26 22:53:22 Mike Exp \$
% Amgen Inc.
% Department of Computational Biology
% mboedigh@amgen.com
%

% todo
%   a simple gui for specifying a model
%   assumption is made that intercept is present in the model, or that it
%   is never explicitly tested
%   does not work as expected with arbitraty model, e.g. [0 0; 1 0; 1 1];

% number of tests
ntests = length(reference);

% get term names from full model
hasIntercept = glm.hasIntercept;
tn = model2eqn(glm.model, glm.var_names, '*');
if hasIntercept
tn(1) = [];
end;

% reference models
rmodel = tests(reference,:);

% test (reduced) models
tmodel = tests(1:ntests,:);

% terms dropped relative to the reference model
dropped_terms = ~rmodel - ~tmodel;

% names of terms left out of test but that are in reference model
dropped_names =  model2eqn( dropped_terms, tn, ',', '' );

% get names of terms in both the test and reference models.
% these are prior terms in that we are testing the effect of the
% dropped terms given that these prior terms are in the model
%
prior_terms   = double(~tmodel&~rmodel);
prior_names = model2eqn( prior_terms, tn, ',' );
prior_names = strcat('|',prior_names);

% get nested terms
% these terms are in dropped in both the reference and test models
nested_terms = tmodel&rmodel;
nested_names = model2eqn( nested_terms, tn, ',', '' );
i = any(nested_terms,2);
nested_names(i) = strcat( '(', nested_names(i), ')');

source = strcat( dropped_names, nested_names, prior_names );

if nargin > 3  && short_form == 1
source = regexprep( source, '\|.*\$', '' );
end;
```