In[]:=
beta:=(Times@@(Gamma[#]))/(Gamma[Total[#]])&;dirichlet:=(Times@@(#1^(#2-1)))/(beta[#2])&;GROUPS:=With[{vars=#},ImplicitRegion[Total@vars==1&&And@@Thread[vars>0],vars]]&
In[]:=
(*setdimensiona=3*)a=3;randVars=Array[x,a];params=θConstantArray[1/a,a];(*createDirichletPDF*)dir=dirichlet[randVars,params];(*comparetotruePDFusingDirichletDistribution*)piecewisePDF=PDF[DirichletDistribution[params],Most@randVars];pdf=piecewisePDF[[1,1,1]];(*subin1-x[1]-...x[a-1]=x[a]*)pdfSubbed=pdf/.((1-Total@Most@randVars)->Last@randVars);(*confirmequality*)pdfSubbed==dir
Out[]=
True
In[]:=
(*getimplicitregionfromDirichletDistributionpdf*)ir1=ImplicitRegion[piecewisePDF[[1,1,2]],#]&@Most@randVars
Out[]=
ImplicitRegion[x[1]>0&&x[2]>0&&1-x[1]-x[2]>0,{x[1],x[2]}]
In[]:=
(*asexpected,weget1whenintegrating*)Integrate[pdf,Most@randVars∈ir1]
Out[]=
1
In[]:=
(*createcustomimplicitregionandintegrate*)ir2=GROUPS@randVars
Out[]=
ImplicitRegion[x[1]+x[2]+x[3]1&&x[1]>0&&x[2]>0&&x[3]>0,{x[1],x[2],x[3]}]
In[]:=
(*wenowget
3
insteadof1*)Integrate[dir,randVars∈ir2]Out[]=
3
In[]:=
(*theissueisthatir2hasaregionmeasurewhiletheconditioninthepiecewisehaspdfhasaregionmeasureof1/2*)RegionMeasure@ir2RegionMeasure@ir1
3
2
Out[]=
3
2
Out[]=
1
2
In[]:=
(*wecanofcoursesubinx[a]->1-Total@Most@randVars,butitsoundslikeyouwouldliketoexplicitlykeepallofthevariablesinthepdf*)subbedDir=dir/.x[a]->1-Total@Most@randVars;Integrate[subbedDir,Most@randVars∈ir1]
Out[]=
1
In[]:=
(*wecanalsoplotrandompointsfromeachregionandseetheyareindeedthesame(thiswilltake10-20secondstorun).Iexplictlysettheseedtogetthesamerandompointseachtimethisisrun,butyoucangetridoforchangetheSeedRandomifyouwanttolookatotherrandompoints*)SeedRandom[1234];n=10^3;(*grabpointsfromregiongeneratedbyDirichletDistribution.*)randPoints1=RandomPoint[ir1,n];(*Ijoin1minusthetotaloftheotherrandomvariablesforthevalueofthelastrandomvariabletomatchthedimensionoftheotherregion*)randPoints1=Join[#,{1-Total[#]}]&/@randPoints1;(*sinceir2isinfinitelythininoneofitsdimensions(since1==Total[randVars]),wehavetomanipulatetheregionslightlyforthealgorithmsusedbyRandomPointtoworkcorrectly,sowerelaxthecondition1==Total[randVars]tobeAbs[1-Total[randVars]]<dforsomesmalld*)relaxedRegion=With[{vars=randVars,d=10^-5},ImplicitRegion[Abs[1-Total[randVars]]<d&&And@@Thread[vars>0],vars]];randPoints2=RandomPoint[relaxedRegion,n];(*thetworegionsappeartobethesame*)ListPointPlot3D[{randPoints1,randPoints2}]
Out[]=
In[]:=
(*ifweonlygrabthefirst2randomvariablesfromeachrandompointlist,theyshouldmatchaswell*){randPoints12D,randPoints22D}=(Most/@#&)/@{randPoints1,randPoints2};ListPlot[{randPoints12D,randPoints22D},PlotRange->All]
Out[]=