1

Let us say I have generated a table tab with entries {p1,..,pd,chi2[p1,..,pd]} of dimension {ngrid[1],..,ngrid[d],d+1} where chi2 is a function of d parameters. More precisely, tab is generated using the following code:

tab = Table[
{p1, p2, p3, p4, p5} = {pTab[1][[i1]], pTab[2][[i2]], pTab[3][[i3]],pTab[4][[i4]], pTab[5][[i5]]};
{p1, p2, p3, p4, p5, chi2[p1, p2, p3, p4, p5]}
, {i1, 1, ngrid[1]}, {i2,1,ngrid[2]}, {i3,1,ngrid[3]}, {i4,1,ngrid[4]}, {i5,1,ngrid[5]}]

where pTab[i] are the tables giving the pi parameters.

Now, I want to marginalize the table over some parameters. More precisely I want to marginalize the corresponding likelihood ($\exp(-\chi^2/2)$). This can be done over all parameters but p2 in the following way:

newtab = tab[[1, All, 1, 1, 1, {2, 6}]];
Do[
newtab[[j,2]]=-2Log[Total[Exp[-tab[[All,j,All,All,All,6]]/2],4]];
,{j, 1,ngrid[2]}];

If one wants to marginalize over all parameters but p1 and p3:

newtab = tab[[All, 1, All, 1, 1, {1, 3, 6}]]
Do[
newtab[[j, k, 3]] = -2Log[Total[Exp[-tab[[j,All,k,All,All,6]]/2],4]];
,{j, 1, ngrid[1]}, {k, 1, ngrid[3]}];

and so on.

It would be great to have a function that can do the tasks above and which works for any number d of parameters and which can accept any combination of parameters to marginalize over. Something like marginalize[{p2,p4,p5},tab] where tab is my table and {p2,p4,p5} are the parameters to marginalize over. Do you know how to code such a function?

Valerio
  • 1,982
  • 13
  • 23
  • I'm not sure whether you're asking about the math or about how to code it in Mathematica – Dr. belisarius Jul 07 '13 at 15:02
  • 1
    Maybe you could provide us with a simple example: say two ro three dimensions and a clear explanation of what the "marginalization" is supposed to do. For example, it is not clear whether this is much more than just normalizing the column/row of the matrix. – bill s Jul 07 '13 at 15:20
  • marginalization is supposed to do exactly what shown in the examples. The problem is that I don't know how to code a function that can do flexibly those jobs. – Valerio Jul 07 '13 at 16:18
  • But Valero, we don't know what your function does because you haven't defined all the relevant data: what is tab? What is ngrid? If you provide an example where it is clear what you are after, you will be more likely to get help. – bill s Jul 07 '13 at 20:41
  • @bill s I added some explanations, is it now clear? – Valerio Jul 18 '13 at 12:31
  • Valerio could you post an example, for me to understand how the program works? Thanks. – user11892 Jan 22 '14 at 01:22

1 Answers1

4

I propose this function

marginalize[pp_, tab_] := 
 Module[{newtab, id, d = Dimensions[tab][[-1]] - 1},
  id = ConstantArray[1, d];
  id[[pp]] = All;
  newtab = Transpose[tab[[Sequence @@ id, {Sequence @@ pp, -1}]], Ordering[pp]];
  With[{nn = Sequence @@ Table[Unique["n"], {Length[pp]}]},
   With[{lim = Sequence @@ Transpose@{{nn}, Dimensions[tab][[pp]]},
     mm = (id = ConstantArray[All, d]; id[[pp]] = {nn}; Sequence @@ id)},
    Do[newtab[[nn, -1]] = -2 Log[Total[Exp[-tab[[mm, -1]]/2], Infinity]];, lim];
    ]];
  newtab
  ];

It does not depend on any external parameters like the number of dimensions, ngrid, etc.

I think tab can be generated by

d = 5;
ClearAll[pTab];
Do[pTab[i] = Range[-1, 1, 0.1], {i, d}];
chi2 = Total[#^2] &; (* I don't know your chi2 and pTab *)

tab = Module[{T, n}, With[{F = chi2}, Replace[
     Hold[Array[Append[T, F[T]] &, Length[pTab[#]] & /@ Range[d]]] /.
        T -> Table[Hold[pTab[n][[Slot[n]]]] /. n -> i, {i, d}],
     Hold[x_] :> x, {0, Infinity}, Heads -> True]]];

tab // Developer`PackedArrayQ

True

It takes the number of dimensions d as parameter. It's quite complicated but it generate packed array, which is much faster.

ybeltukov
  • 43,673
  • 5
  • 108
  • 212
  • thanks! only one thing: the parameters pp have to be sorted. Or, would it be possible to modify your function such that if i give pp={3,2}, then i get a table which has the third parameter in the first dimension and the second parameter in the second dimension? Otherwise i can do this later using: newtab = newtab[[All, All, {2, 1, 3}]]; newtab = Transpose[newtab, {2, 1}] – Valerio Sep 12 '13 at 10:01
  • @Valerio I add transposition to my code. I think it now supports unsorted pp. Please check it. – ybeltukov Sep 12 '13 at 16:09
  • It's working alright, thanks! – Valerio Sep 13 '13 at 10:15
  • could you please parallelize marginalize? Perhaps Do -> ParallelDo? – Valerio Apr 04 '18 at 18:55