function [] = sumofsquares(d)

n = elems(d);

ssy = ssq(d,n)
gm = matmean(d);

ss0 = gm*gm*n
sst = ssy - ss0

% sum of squares main effects
e = effects(d);
e2 = e .* e;
ssa = (n/3) * sum(e(1,:) .* e(1,:))
ssd = (n/3) * sum(e(2,:) .* e(2,:))
ssp = (n/3) * sum(e(3,:) .* e(3,:))
ssm = (n/3) * sum(e(4,:) .* e(4,:))

% sum of squares for first-order interactions

[ad, ap, am, dp, dm, pm] = interact1(d);
dm
ssad = ssq(ad,n)
ssap = ssq(ap,n)
ssam = ssq(ap,n)
ssdp = ssq(dp,n)
ssdm = ssq(dm,n)
sspm = ssq(pm,n)
firstorderssq = (ssad + ssap + ssam + ssdp + ssdm + sspm)
firstorderpv = 100 *firstorderssq/sst

% sum of squares for second-order interactions

[adp, adm, apm, dpm] = interact2(d);
ssadp = ssq(adp,n)
ssadm = ssq(adm,n)
ssapm = ssq(apm,n)
ssdpm = ssq(dpm,n)
secondorderssq = (ssadp + ssadm + ssapm + ssdpm)
secondorderpv = 100 *secondorderssq/sst

% sum of squares for third-order interactions
[adpm] = interact3(d);
ssadpm = ssq(adpm,n)
thirdorderssq = (ssadpm)
thirdorderpv = 100 *thirdorderssq/sst
