# Dataset and Select with a SmoothKernelDistribution

This question is hopefully a case of stupidity on my side, but I just don't understand why this happens. I had code that worked perfectly fine for several weeks and gives an error now on Linux and Mac. I could reduce it to a simple example.

First, we create a simple Dataset consisting of 2d random points and I create a kernel density estimator for the point set:

RandomSeed[12345];
ds = Dataset[<|"a" -> #1, "b" -> #2|> & @@@
RandomVariate[
MultinormalDistribution[{24.6, 104.19}, {{10.2, 1}, {1, 15.9}}],
407]];
kd = SmoothKernelDistribution[Normal[ds[All, Values]]]


Now I want to select some of the points that are inside the area of a certain probability which can be done with the PDF and a simple Select:

kernelProbability = 0.00059914;
func = PDF[kd, {#"a", #"b"}] >= kernelProbability &;
ds[Select[func] /* Length]


What I get is an unspecific error message Part::partd

Turning the call inside out and using Select[ds, func] doesn't change a thing. I was debugging this by trying whether my func did work, and it does

func@Normal[ds[1]]
(* True *)


Here comes the crazy part, if you now evaluate the Select call again, it does work

ds[Select[func] /* Length]
(* 397 *)


I think I could trace the source of the error down, but I don't see why this happens. Does anyone have a clue? Can you reproduce this behavior?

# Update

I reported the issue and it is now filed for the Wolfram dev team: CASE:3826715

• I can reproduce this behaviour in V11.0.1 on Win7/64. ds[1 /* func] fails but ds[1 /* func, FailureAction->None] succeeds and causes the first expression to start working. Neither ClearSystemCache[] nor DatasetResetCompileQueryCache[] returns the system to the broken state. If I define vals = Normal[ds[All, Values]] then evaluating kd = SmoothKernelDistribution[vals] alone is enough to restore the broken behaviour. Jan 25, 2017 at 15:47
• I can reproduce on 10.4.1 as well. Jan 25, 2017 at 16:32
• Thank god I'm not that crazy. I may have figured out why my code worked the whole time. The example I posted is only a small part of a larger notebook. Yesterday, I was working only on that particular section without evaluating the whole notebook. It seems that at some point, I did something like func@Normal[ds[1]] to test stuff and therefore, the issue never appeared before. Very weird. Jan 26, 2017 at 4:06

Short, Short Version

A General::stop message will cause a dataset query to fail if some other explicitly quieted message occurs three or more times. For example:

Dataset[1][Quiet[Do[Message[f::ivar, #], 2]; #, f::ivar] &]
(* 1                -- success: only two quieted f::ivar messages generated *)

Dataset[1][Quiet[Do[Message[f::ivar, #], 3]; #, f::ivar] &]
(* Failure[...]     -- third f:ivar message triggers a General::stop and failure *)


Short Version

The behaviour we see is due to an unfortunate combination of circumstances:

1. Generating the PDF of a SmoothKernelDistribution generates some messages which are internally suppressed using Quiet.
2. Distribution results are routinely cached so that the messages in question are only generated during the first evaluation.
3. FailureAction -> "Abort" is the default setting for Dataset queries and it causes queries to fail if any unexpected messages are generated. A message is unexpected if it has not been explicitly quieted by code executed within the body of the query.
4. PDF[...] is generating so many Part::partd messages that the system issues a General::stop message. Normally, a General::stop receives special treatment and quieted along with the messages that generated it, but the hypersensitive FailureAction machinery classifies it as an unexpected message and fails the query.

Work-around

A work-around is to turn off the special message handling from the Dataset query:

ds[Select[func] /* Length, FailureAction -> None]
(* 397 *)


Is It A Bug?

Yes, I think so. Given that in normal circumstances General::stop is only generated after some other message has appeared for the third time, the bug would be corrected if General::stop was always presumed to be quieted by GeneralUtilitiesMessageQuietedQ (see discussion below).

Spelunking Report (current as of version 11.0.1)

... be prepared for too much complicated trivia ...

Let's begin by finding a minimal set of steps to reproduce the behaviour:

d = SmoothKernelDistribution[{{1.,2.}}];

Dataset[{1., 2.}][PDF[d, #]&]
(* Failure[...] *)

PDF[d, {1., 2.}]
(* 0.390901 *)

Dataset[{1., 2.}][PDF[d, #]&]
(* 0.390901 *)

d = SmoothKernelDistribution[{{1.,2.}}];

Dataset[{1., 2.}][PDF[d, #]&]
(* Failure[...] *)


Notice how the bad behaviour disappears after performing a non-query evaluation and how it reappears after reassigning a newly created distribution to d. (bonus trivia: using the operator form PDF[d] always fails as it appears to defeat the cache)

PDF Evaluation Generates Messages

By blocking the action of Quiet, we can observe that PDF[d, {1., 2.}] generates messages that would normally be quieted:

d = SmoothKernelDistribution[{{1.,2.}}];
Block[{Quiet}, PDF[d, {1., 2.}]]

(*
>> Part::partd: Part specification $x[[1]] is longer than depth of object. >> Part::partd: Part specification$x[[2]] is longer than depth of object.
>> Part::partd: Part specification \$x[[1]] is longer than depth of object.
>> General::stop: Further output of Part::partd will be suppressed during this calculation.

0.390901
*)


... but only the first time:

Block[{Quiet}, PDF[d, {1., 2.}]]
(* 0.390901 *)


We can use an alternative technique to view considerably more details about the generated messages and their evaluation contexts:

d = SmoothKernelDistribution[{{1.,2.}}];
InternalHandlerBlock[{"Message",Print[InternalQuietStatus[]]&}, PDF[d,  {1., 2.}]]

(*
{Global->Unquiet,Off->{Part::partd},On->{},Stack->{{1458,{Part::partd},{}}},MessageList->{{Part::partd,1458}},Check->None}
{Global->Unquiet,Off->{Part::partd},On->{},Stack->{{1458,{Part::partd},{}}},MessageList->{{Part::partd,1458},{Part::partd,1458}},Check->None}
...
{Global->Unquiet,Off->{Part::partd},On->{},Stack->{{1458,{Part::partd},{}}},MessageList->{{Part::partd,1458},{Part::partd,1458},{Part::partd,1458},{General::stop,1458}},Check->None}
{Global->Unquiet,Off->{Part::partd},On->{},Stack->{{1458,{Part::partd},{}}},MessageList->{{Part::partd,1458},{Part::partd,1458},{Part::partd,1458},{General::stop,1458}},Check->None}

0.390901
*)


We will shortly see that the Stack values from InternalQuietStatus are relevant.

Distribution Calculations Use A Cache

The presence of a cache is strongly hinted by the fact that the messages disappear during the first evaluation of PDF but reappear if we regenerate the distribution itself. We can confirm this guess by tracing calls to StoreDataDistributionExpression:

On[StatisticsDataDistributionUtilitiesStoreDataDistributionExpression]

SystemDumpDeactivateReadProtected[{DataDistribution, PDF}
, Print["#### define"]; d = SmoothKernelDistribution[{{1.,2.}}]
; Print["#### pdf1"];   PDF[d,{1.,2.}]
; Print["#### pdf2"];   PDF[d,{1.,2.}]
]

Off[]

(* Output:
#### define
... messages showing calls to StoreDataDistributionExpression ...
#### pdf1
... messages showing more calls to StoreDataDistributionExpression ...
#### pdf2
... no messages! ...
*)


Observe how the initial distribution generation tucked some values into a cache, as did the first PDF evaluation. But the second PDF evaluation did not store into the cache at all.

FailureAction Treats Messages As Query Failure

Now, we can move on to FailureAction in dataset queries. As noted earlier, the default action is to cause a query to fail should any messages be generated:

Dataset[1][(Message[f::ivar, #]; #) &]
(* Failure[...] *)


But quieted messages do not (normally) cause failure:

Dataset[1][Quiet[Message[f::ivar, #]; #] &]
(* 1 *)


We can also turn off the failure processing:

Dataset[1][(Message[f::ivar, #]; #) &, FailureAction -> None]
(*
>> f: 1 is not a valid variable
1
*)


Unexpected Messages Trigger FailureAction

But if the previous example shows that quieted messages do not trigger FailureAction, what is so special about the quieted messages we observed for this:

d = SmoothKernelDistribution[{{1.,2.}}];
Dataset[{1., 2.}][PDF[d, #]&]
(* Failure[...] *)


To find the answer, we need to closely inspect a voluminous trace of the evaluation (not reproduced here). Deep in the bowels of that trace, we find that the relevant functions are:

Needs["GeneralUtilities"]
{ EvaluateChecked, GeneralUtilitiesFailurePackagePrivateCheckedHandler
, MessageStackID, MessageQuietedQ
} // Scan[PrintDefinitionsLocal]


In short, EvaluateChecked begins by setting up an environment which will catch and handle any messages. MessageStackID is used to identify the outermost stack boundary of that environment. Should a message appear, CheckedHandler will inspect the message to determine if it should be ignored. Any message that is not MessageQuietedQ will cause a failure. MessageQuietedQ uses the InternalQuietStatus[] output we saw earlier to see if the message has been explicitly quieted by code within the bounded stack environment.

Putting It All Together

In the case at hand, PDF[...] explicitly quiets the message Part::partd. But is does not explicitly quiet General::stop.

And so, to finally reach the end of our long-winded shaggy dog story... it is the General::stop message that causes the query to fail. That message only appears the first time the PDF function is generated because the function is cached on subsequent attempts.

• And it explains perfectly, why I haven't noticed this earlier. I was of course inspecting the PDF and creating some density plots where I evaluated the PDF. Thanks for this answer. Jan 28, 2017 at 5:41
• Hell of a shaggy dog indeed. +1! Jan 28, 2017 at 6:13
• @halirutan My pleasure -- I was intrigued by this strangeness of this one. Thanks for the bounty. Jan 30, 2017 at 15:18

# Workaround with Efficiency Improvement

Another approach, which also makes your code more efficient, is to create the PDF of the SmoothKernelDistribution only once.

When written as

func = PDF[kd, {#"a", #"b"}] >= kernelProbability &;
ds[Select[func] /* Length]


kd is resolved into its PDF for each row in the Dataset because, as Attributes tells us, Function has attribute HoldAll. Here it is not too complicated for PDF to resolve kd into its PDF so you don't really notice this.

However, we can resolve the PDF of kd only once and then evaluate it repeatedly in the Query.

ClearAll[x, y];
pdf[x_, y_] = PDF[kd, {x, y}];
func1 = pdf[#"a", #"b"] >= kernelProbability &;
ds[Select[func1] /* Length]


pdf is Set to the resolved PDF of kd outside of the Select. Now for each row the resolved PDF (pdf) is evaluated. The overhead of resolving it for each row is removed. This not only makes your code more efficient but bypasses the bug as well.

In general be aware of the cost of your expression when doing row-by-row evaluations.

Hope this helps.