Here is a specific example of how to do it. If you need to generalize it and can't figure that out, let me know. The idea is to multiply the first column of d1 by each column of d2, then the second column of d1 by each column of d2:
d1 = dummyvar(T); d2 = dummyvar(R);
d1(:,[1 1 1 2 2 2]).*d2(:,[1 2 3 1 2 3])
Another thing to consider is the x2fx function, using something like this:
x2fx([T R],'interaction',1:2)
This requests that both T and R be treated as categorical, and that all columns for an interaction model be computed. You'll see the result has a column of ones, then dummy columns for T and R separately, then columns for the interaction. The idea is that these will be used in regression, so certain columns are omitted so that the full matrix is not overdetermined.
Best Answer