Advertisement
Not a member of Pastebin yet?
Sign Up,
it unlocks many cool features!
- pdf = HistogramDistribution[data];
- DensityPlot[PDF[pdf, {x, y}], {x, 1.5, 3.3}, {y, -46, -44}, PlotPoints -> 100]
- base = RandomVariate[BinormalDistribution[{2.4, -44.9}, {0.25, 0.3}, 2/3], 400];
- contam = RandomVariate[BinormalDistribution[{2.4, -45.3}, {0.1, 0.2}, -0.9], 100];
- data = base~Join~contam;
- ListPlot[{base, contam}]
- pdf = SmoothKernelDistribution[data, 0.1, "Gaussian"];
- plot = ContourPlot[PDF[pdf, {x, y}], {x, 1.5, 3.3}, {y, -46, -44},
- PlotRange -> {Full, Full, Full}, Contours -> 23, ColorFunction -> "DarkRainbow"]
- f[x_, y_, t_: 0] := With[{z = PDF[pdf, {x, y}]}, z Boole[z >= t]];
- r = FindRoot[NIntegrate[f[x, y, t], {x, 1.5, 3.3}, {y, -46, -44},
- AccuracyGoal -> 3, PrecisionGoal -> 6] - 0.9, {t, 0, 2},
- Method -> "Brent", AccuracyGoal -> 3, PrecisionGoal -> 6]
- NIntegrate[f[x, y, t /. r], {x, 1.5, 3.3}, {y, -46, -44}] /. r
- ci = ContourPlot[PDF[pdf, {x, y}] == (t /. r), {x, 1.5, 3.3}, {y, -46, -44},
- ContourStyle -> Directive[Thick, White]];
- Show[plot, ci]
- Needs["MultivariateStatistics`"];
- data=RandomVariate[BinormalDistribution[{2.4,-45},{.34,.31},.3],10^4];
- dist=HistogramDistribution[data,"FreedmanDiaconis"];
- conf=EllipsoidQuantile[data,.9];
- confConv=PolytopeQuantile[data,.9];
- Show[DensityPlot[Evaluate@PDF[dist, {x, y}], {x, 1.5, 3.3}, {y, -46, -44},
- PlotPoints -> 100, ColorFunction -> "Pastel", Exclusions -> None,
- PlotRange -> All],
- Graphics[{Directive[LightOrange, Opacity[.4]],
- EdgeForm[Directive[Orange, Opacity[.6], Thick]],
- FilledCurve[BSplineCurve[Graphics[confConv][[1, 1]], SplineClosed -> True]]}],
- Graphics[{Directive[Black, Opacity[.7]], Thick, Dashed, conf}],PlotRange -> All]
- weight = RandomVariate[NormalDistribution[10, 2.43], 10^4];
- weightedData = WeightedData[data, weight];
- smdistwt = HistogramDistribution[weightedData, "FreedmanDiaconis"];
- dataWeighted = RandomVariate[smdistwt, 10^4];
- confwt = EllipsoidQuantile[dataWeighted, .9];
- confConvwt = PolytopeQuantile[dataWeighted, .9];
- Needs["Combinatorica`"] (* for BinarySearch *)
- findCredibilityLevel[data_, level_] := Module[{srtdata, cumsum, index},
- srtdata = data // Flatten // Sort;
- cumsum = Accumulate[srtdata];
- (*overcover interval with Floor, Ceil would undercover*)
- index = BinarySearch[cumsum, (1 - level)*Last[cumsum]] // Floor;
- srtdata[[index]]
- ]
- xmin = 1.5; xmax = 3.3; ymin = -46; ymax = -44; npix = 100;
- base =
- RandomVariate[BinormalDistribution[{2.4, -44.9}, {0.25, 0.3}, 2/3],
- 400];
- contam = RandomVariate[
- BinormalDistribution[{2.4, -45.3}, {0.1, 0.2}, -0.9], 100];
- data = base~Join~contam;
- pdf = SmoothKernelDistribution[data, 0.1, "Gaussian"];
- plot = ContourPlot[PDF[pdf, {x, y}], {x, xmin, xmax}, {y, ymin, ymax},
- PlotRange -> {Full, Full, Full}, Contours -> 23,
- ColorFunction -> "DarkRainbow"];
- dx = (xmax - xmin)/npix;
- dy = (ymax - ymin)/npix;
- data = Table[
- PDF[pdf, {xmin + i dx, ymin + j dy}], {i, 1, npix}, {j, 1, npix}];
- t = findCredibilityLevel[data, 0.9];
- ci = ContourPlot[
- PDF[pdf, {x, y}] == t, {x, xmin, xmax}, {y, ymin, ymax},
- ContourStyle -> Directive[Thick, White]];
- Show[plot, ci]
Advertisement
Add Comment
Please, Sign In to add comment
Advertisement