Moment vs. loglinear polytope

Init

walshFeatures:=Module[{},
Times@@@(Table[Subsets[vars,{order}],{order,0,d}]//Flatten[#,1]&)
];

(* Give design matrix for given feature functions *)
designMatrix:=(
Table[feats[[i]]/.Thread[vars->xs[[j]]],{i,Length[feats]},{j,Length[xs]}]
);

masks[len_,num_]:=Module[{},
arr=Array[0&,len];
Table[ReplacePart[arr,#->1&/@t],{t,Subsets[Range[len],{num}]}]
]
initDomain[numvars_]:=(
d=numvars;
dd=(2^d);
vars=Array[xx,d];
feats=walshFeatures;

domain={1,-1};
xs=Tuples[domain,d];
DD=designMatrix;
EE=Inverse[designMatrix];
p=Array[pp,dd];

u=Array[1&,dd];

P$p=p/Total[p];
P$μ=EE.p;
P$θ=Exp[p.DD]/Total[Exp[p.DD]];

(* three parameterizations. i=1 is p, i=2 is mu, i=3 is theta.
conv[[i1,i2]] gives values of i1 in terms of i2 *)
(* ie, conv[[2,1]] will convert p into mu *)
conv=(({{p, EE.p+1/dd, Exp[p.DD]}, {DD.(p-1/dd), p, DD.Exp[p.DD]}, {Log[p].EE, Log[EE.p].EE, p}}));

(* how to show each parameterization on screen *)
(*screen={Rest[DD.p],Rest[p],Rest[p]};*)
screen={Rest[DD.(p/Total[p])],Rest[p],Rest[p]};
center={u/dd,u*0,u*0};
range={{0,3},{0,1},{0,2}};
ev[expr_,vals_]:=expr/.Thread[p->vals];
);

initDomain[2];
simplex=Graphics3D[{Yellow,Opacity[.2],GraphicsComplex[Transpose@Rest@DD,Polygon[Subsets[{1,2,3,4},{3}]]]}];

Use

(* generates types for dice with dims dimensions, n samples, min samples in each category *)
getTypes[n_,dims_,min_]:=
1/n Permutations/@IntegerPartitions[n,{dims},Range[min,n]]//Flatten[#,1]&;


types=getTypes[20,4,1];
mu=Rest[ev[conv[[2,1]],#]]&/@types;
theta=Rest[ev[conv[[3,1]],#]]&/@types;

gr1=List@@(ComputationalGeometry`Methods`ConvexHull3D[theta,Graphics`Mesh`FlatFaces->False]);

regular=Show[simplex,Graphics3D[{PointSize[.01],Point[mu]}]];
logSimplex=Graphics3D[{FaceForm[Opacity[.5]],EdgeForm[Opacity[0]]}~Join~gr1];
nonregular=Show[logSimplex,Graphics3D[{PointSize[.005], Point[theta]}]];
GraphicsRow[{regular,nonregular},ImageSize->600]