Comment créer un graphique de distribution marginale en utilisant DensityPlot?

Sep 09 2020

Comment feriez-vous un diagramme de distribution marginale en utilisant DensityPlot? En adaptant la réponse de Sjoerd à cette question en utilisant Epilog, je peux tracer une tranche à une coordonnée donnée (par exemple, x = 0 ou y = 0) comme ci-dessous:

f[x_, y_] := Exp[-2 (x^2 + y^2)] HermiteH[2, Sqrt[2] x]^2;
DensityPlot[f[x, y], {x, -2, 2}, {y, -2, 2}, PlotRange -> All, Frame -> False,
    Epilog -> {Line[Table[{x1, 2.05 + 0.2 f[x1, 0]}, {x1, -2, 2, 0.01}]],
              Line[Table[{ 2.05 + 0.2 f[0, y1], y1}, {y1, -2, 2, 0.01}]] },
    PlotRangePadding -> 0, PlotRangeClipping -> False, ImagePadding -> {{0, 100}, {0, 100}}]

Cela donne le résultat suivant:

Mais ce que j'aimerais vraiment, c'est tracer les valeurs intégrées de colonne (ligne) le DensityPlotlong des marges de l'axe x (y).

La fonction réelle d'intérêt pour ce calcul est coûteuse, il serait donc préférable de n'évaluer qu'une seule fois (c'est-à-dire, pouvons-nous utiliser les valeurs de DensityPlot?). De plus, pour cette raison, il DensityPlotest préféré ListDensityPlotà son échantillonnage automatique de maillage, car les fonctions d'intérêt ont tendance à être localisées, de sorte qu'un maillage uniforme serait inutile.

Voici un exemple où le traçage le long d'une tranche donnée n'est pas représentatif de la distribution marginale:

et les graphiques marginaux (en rouge et bleu, pris le long des lignes pointillées claires) évaluent simplement à zéro.

Merci d'avance!!

Réponses

3 kglr Sep 10 2020 at 05:27

ClearAll[f, xMargin, yMargin, ppX, ppY]

f[x_, y_] := Exp[-2 (x^2 + y^2)] HermiteH[2, Sqrt[2] x]^2

xMargin[x_] = Integrate[f[x, y], {y, -Infinity, Infinity}];

yMargin[y_] = Integrate[f[x, y], {x, -Infinity, Infinity}];

xrange = {-3, 3};

yrange = {-2, 2};

scale = 1/4/Pi;

gap = 0.05;

dp = DensityPlot[f[x, y], {x, xrange[[1]], xrange[[2]]}, {y, yrange[[1]], yrange[[2]]}, 
       PlotRange -> All]

Nous pouvons construire des marges correctement traduites en utilisant ParametricPlot:

ppY = ParametricPlot[{xrange[[1]] - gap - scale v yMargin[y], y},
     {y, yrange[[1]], yrange[[2]]}, {v, 0, 1},
   PlotStyle -> Red, PlotPoints -> 50, Axes -> False];

ppX = ParametricPlot[{x, yrange[[1]] - gap - scale v xMargin[x] },
   {x, xrange[[1]], xrange[[2]]}, {v, 0, 1}, 
 PlotStyle -> Blue, PlotPoints -> 50, Axes -> False];

et combinez-les en dputilisant Show:

Show[ppY, ppX, dp, PlotRange -> All, Frame -> True]

Pour afficher les tracés marginaux sur les cadres supérieur et droit:

ppY2 = ParametricPlot[{xrange[[2]] + gap + scale v yMargin[y], y}, 
    {y, yrange[[1]], yrange[[2]]}, {v, 0, 1}, 
    PlotStyle -> Red, PlotPoints -> 50, Axes -> False];

ppX2 = ParametricPlot[{x, yrange[[2]] + gap + scale v xMargin[x]}, 
    {x, xrange[[1]], xrange[[2]]}, {v, 0, 1}, 
    PlotStyle -> Blue, PlotPoints -> 50, Axes -> False];

Show[ppY2, ppX2, dp, PlotRange -> All, Frame -> True]

Pour mettre les parcelles marginales en dehors du cadre, nous pouvons utiliser Inset+ Epilog:

insetY = Inset[#, {xrange[[2]] (1 + gap), yrange[[2]]},
   {Left, Top}, Scaled[1]] & @ ppY2;

insetX = Inset[#, {xrange[[2]], yrange[[2]] (1 + gap)}, 
   {Right, Bottom}, Scaled[1]] & @ ppX2;

Show[dp, Epilog -> {insetX, insetY}, 
 ImagePadding -> {{Scaled[.02], Scaled[.1]}, {Scaled[.02], Scaled[.1]}}, 
 ImageSize -> Large, PlotRangeClipping -> False, ]

Alternativement, nous pouvons Plotles fonctions xMarginet yMarginet utiliser GeometricTransformationavec les fonctions de transformation appropriées les positionner et Showles objets graphiques transformés avec dp:

ClearAll[transform, tFX, tFY]
transform[tf_] := Graphics[#[[1]] /. 
     ll : (_Line | _Polygon) :> GeometricTransformation[ll, tf]] &;

tFY = TranslationTransform[{-gap, xrange[[1]]}]@*
   RotationTransform[Pi/2, {xrange[[1]], 0}];

tFX = TranslationTransform[{0, yrange[[1]] - gap}]@*
   ScalingTransform[{1, -1}];

pltY = Plot[scale yMargin[y], {y, yrange[[1]], yrange[[2]]}, 
   Filling -> Axis, PlotStyle -> Red, Axes -> False];
pltX = Plot[scale  xMargin[x], {x, xrange[[1]], xrange[[2]]}, 
   Filling -> Axis, PlotStyle -> Blue, Axes -> False];

Show[transform[tFY]@pltY, transform[tFX]@pltX, dp, PlotRange -> All, 
 Frame -> True]

Pour afficher les tracés marginaux sur les cadres supérieur et droit, utilisez les transformations tFX2et tFY2:

tFY2 = TranslationTransform[{gap, xrange[[1]]}]@*
   RotationTransform[-Pi/2, {xrange[[2]], 0}];

tFX2 = TranslationTransform[{0, yrange[[2]] + gap}];

Show[transform[tFY2] @ pltY, transform[tFX2] @ pltX, dp, PlotRange -> All,
  Frame -> True]

Mise à jour: Une approche alternative pour obtenir les tracés marginaux: Utilisez Plot3Dpour tracer favec des lignes de maillage également espacées dans les directions x et y et extraire les coordonnées des lignes de maillage.

ndivs = 50;
{meshx, meshy} = Subdivide[#[[1]], #[[2]], ndivs] & /@ {xrange, yrange};

coords = Plot3D[f[x, y],
    {x, xrange[[1]], xrange[[2]]}, {y, yrange[[1]],  yrange[[2]]}, 
  PlotRange -> All, Mesh -> {meshx, meshy}, PlotStyle -> None][[1, 1]];

Regroupez coordspar les première et deuxième coordonnées et construisez deux WeightedDataobjets et tracez-les en utilisant SnoothHistogram:

bw = .01;
{wDx, wDy} = Table[Apply[WeightedData] @ Transpose @ KeyValueMap[List] @
      KeySort @ GroupBy[coords, Round[#[[i]], bw] & -> Last, Mean], {i, 2}];

{sHx, sHy} = {SmoothHistogram[wDx, PlotStyle -> Blue, 
    Filling -> Axis, ImageSize -> 300], 
   SmoothHistogram[wDy, PlotStyle -> Red, Filling -> Axis, ImageSize -> 300]};

Row[{sHx, sHy}, Spacer[10]]

Vous pouvez également Plotl' PDFde SmoothKernelDistributionde wDxet wDy:

{sKDx, sKDy} = SmoothKernelDistribution /@ {wDx, wDy};

{sHx2, sHy2} = {Plot[PDF[sKDx]@x, {x, xrange[[1]], xrange[[2]]}, 
    PlotStyle -> Blue, Filling -> Axis, ImageSize -> 300], 
   Plot[PDF[sKDy]@y, {y, xrange[[1]], yrange[[2]]}, PlotStyle -> Red, 
    Filling -> Axis, ImageSize -> 300]};

Row[{sHx2, sHy2}, Spacer[10]]

Mise à jour 2: Traitement de la DensityPlotsortie pour obtenir des {x,y,z}coordonnées (où zest mis à l'échelle à l'intervalle unitaire: