# Smoothing ListContourPlot contours

I have some data of type {{x1, y1, f1}, {x2, y2, f2}, ...} and want to do ListContourPlot. However, there is a problem that the data f is not smooth enough and has small errors.

Here is a toy example illustrating the problem:

data = Flatten[
Table[{i (1 + RandomReal[0.1]), j (1 + RandomReal[0.1]),
i + j}, {i, 0, 100}, {j, 0, 100}], 1];
ListContourPlot[data, Contours -> {80, 120}, ContourShading -> None]


I would need to smooth the lines in the above curve, to get something like

Any ideas to get this kind of smoothed plot? Thanks!

And here is the figure from real data (which is too large to paste here), which I need to smooth:

I find this thread InterpolationOrder for ContourPlot related. But I was unable to get the method work in my case because of (I guess) different input.

Thanks a lot for the answers! Those answers works great for the toy example, but for the real data I still cannot get smooth lines so far (edit: so far means before the great guys update their answers ^_^. Now it worked great). Here is my data and plotting function, just in case you would like to give a try. I will also investigate why those interesting methods fail to work when I apply it to data...

https://www.dropbox.com/s/2pgobgrhbo42f4v/contourData.dat

data = << "contourData.dat";
ListContourPlot[data, Contours -> {2.30, 6.18},

• I guess the main issue is that your toy example is uniformly sampled, while your real data is very far from it (see ListPlot[Most /@ data, AspectRatio -> 1]: i.stack.imgur.com/DgAAh.png vs. i.stack.imgur.com/zfE5C.png). Also, the aspect ratio of the data is very different from the way you're plotting it. – user484 Apr 7 '14 at 17:21
• @RahulNarain : Thanks! Yes it's a very good point. The real data is from Markov chain Monte Carlo (MCMC), which is by purpose not uniformly distributed. – Yi Wang Apr 7 '14 at 19:43

N.B. Your actual data calls for a more sophisticated approach than the quick hack in my original answer, so I've replaced it with a much better and quite general solution.

There are two things that make your actual data harder to work with than the toy example. First, it is highly irregular and nonuniformly distributed:

ListPlot[Most /@ data, AspectRatio -> 1]


and second, it has an awful aspect ratio:

ListPlot[Most /@ data, AspectRatio -> Automatic]


If your $x$ and $y$ axes aren't actually spatial coordinates but represent independent quantities with different units, you would do well to rescale the data so that, say, the variances along both axes are equal:

sd = StandardDeviation[data];
{scalex, scaley} = 1/Most[sd];
scaledData = {scalex, scaley, 1} # & /@ data;


Now to smooth an arbitrary nonlinear function described by noisy samples at irregularly scattered locations, I believe a local regression (LOESS) model is appropriate. In fact, LOESS is a generally useful thing so it's worthwhile to have an implementation available in Mathematica, and since it's pretty simple to implement, I went ahead and did it. My implementation follows Cleveland and Devlin's 1988 paper, i.e. local quadratic regression with tricube weights, except their $q$ is my $k$.

(* precompute spatial data structure because we'll be needing nearest neighbours a lot *)
nearest = Nearest[scaledData /. {x_, y_, z_} :> ({x, y} -> {x, y, z})];
(* local quadratic regression with k neighbours around point (x, y) *)
loess[nearest_, k_][x_, y_] := Module[{n, d, w, f},
n = nearest[{x, y}, k];
d = EuclideanDistance[{x, y}, Most[#]] & /@ n;
d = d/Max[d];
w = (1 - d^3)^3;
f = LinearModelFit[n, {u, v, u^2, v^2, u v}, {u, v}, Weights -> w];
f[x, y]]


Now we can plot the contours of the regression function, scaling the coordinates back to the original data:

fit = loess[nearest, 100];
{{xmin, xmax}, {ymin, ymax}, {zmin, zmax}} = {Min[#], Max[#]} & /@ Transpose[data];
ContourPlot[fit[scalex x, scaley y], {x, xmin, xmax}, {y, ymin, ymax},
Contours -> {2.30, 6.18}, ContourShading -> None]


Seems to work pretty well.

• Isn't this sort of MovingAverage? – Vitaliy Kaurov Apr 7 '14 at 16:33
• @RahulNarain I would suggest to add smouthData[n_] and this gif for fanciness purposes ;o) – Öskå Apr 7 '14 at 16:43
• The result is amazing! I will read it carefully and try to understand the algorithm! – Yi Wang Apr 8 '14 at 7:24

This is likely a consequence of taking the toy example too seriously, but LinearModelFit seems like a good choice:

lmf=LinearModelFit[data,{x,y},{x,y}]



For a the data provided you might get some use from:

basis[n_,m_] := Flatten[Table[x^i y^j,{i,0,n},{j,0,m}],1]

lmf = LinearModelFit[data,basis[4,6],{x,y}];


Now use this linear model:

cp = ContourPlot[lmf[x, y], {x, 0, 0.5}, {y, -4, 4},
Contours -> {2.30, 6.18}, ContourShading -> None,
ContourStyle -> Red]

Show[ListContourPlot[data, Contours -> {2.30, 6.18},


Not too bad, of course you can adjust the basis used as appropriate. At some point, if you have the information to create a nonlinear model that will be better of course.

• If I understand correctly, his real data are not linear and the model probably is unknown. – Vitaliy Kaurov Apr 7 '14 at 16:31
• NonlinearModelFit might also be appropriate as well. In either case, I suspect they don't want to "smooth" as much as fit. The picture they show doesn't look to terrible to model. – chuy Apr 7 '14 at 16:34
• Sorry, I'm facing a similar problem. How can I use NonlinearModelFit in this algorithm? – Orion Jun 23 at 8:32

I was thinking - how can we average but without loss of points? Well we can randomly sample, interpolate and average - as many times as we want.

Let's take a look at more complicated data:

data = Flatten[Table[{i (1 + RandomReal[0.1]), j (1 + RandomReal[0.1]),
i^2 + j^2}, {i, 0, 100}, {j, 0, 100}], 1];


This is like ~ 10^4 points. Grab samples by 1000 - and many of those - and Interpolate - ListContourPlot anyways does that:

Do[f[k] = Interpolation[RandomSample[data, 1000], InterpolationOrder -> 1], {k, 100}]


Average them:

F[x_, y_] = Sum[f[k][x, y], {k, 100}]/100;


Now let's see:

Show[
ContourPlot[F[x, y], {x, 0, 100}, {y, 0, 100},
ContourShading -> None, ContourStyle -> Directive[Red, Thick]],
ListContourPlot[data, ContourShading -> None, Contours -> 9]
] // Quiet


Anyway - something along these lines.

============ OLD VERSION ============

Well, what I will suggest is brutally simple. You want something better. But for the sake of thought entertainment... You get too much detailed interpolation because you got so many data points. Reduce them?

Manipulate[
Show[
ListPlot3D[data[[1 ;; -1 ;; n]], PlotStyle -> Opacity[.5],
MeshFunctions -> {#3 &}, Mesh -> 10,
MeshStyle -> Directive[Opacity[.5], Thick],
PerformanceGoal -> "Quality"],
ListPointPlot3D[data[[1 ;; -1 ;; n]],
PlotStyle -> Directive[Opacity[.5], Red]]
, ImageSize -> 500]
, {{n, 1}, 1, 147, 1, Appearance -> "Labeled"}]


ListContourPlot[data[[1 ;; -1 ;; 139]], Contours -> {80, 120}, ContourShading -> None]


You can also use moving average, but it also removes points:

ListContourPlot[MovingAverage[data, 50], ContourShading -> None, Contours -> {80, 120}]

• Hey Vitaly, you really are into the screencast thing (I like it) ;D – Yves Klett Apr 7 '14 at 16:22
• @YvesKlett haha yes I love that ;-) – Vitaliy Kaurov Apr 7 '14 at 16:26
• @VitaliyKaurov : Thanks a lot for the answer! The data[[a ;; b ;; c] method indeed improves the figure. But as you said, some data are lost and in my real example the quality is still not good enough. The moving average way is lovely. However, my {{x1, y1, f1}, {x2, y2, f2}, ...} list is actually not ordered in x1 < x2 < x3... but instead in a random way. Thus MovingAverage does not directly apply for my real need... Also, the screencasting is amazing :) – Yi Wang Apr 7 '14 at 16:33
• @YiWang Then RahulNarain moving average method via nearest should work. – Vitaliy Kaurov Apr 7 '14 at 16:35
• @VitaliyKaurov : Thanks a lot for the update! Yes, random sampling, re-sampling and averaging is a very interesting idea :) – Yi Wang Apr 8 '14 at 11:12