Prev: Filter and maximal graylevel
Next: [speed up] vectorized find() for sorted data?
From: Ryan on 2 Feb 2010 12:25
I've got a couple questions about the implementation of loess/lowess in the Matlab smooth function. Specifically, the tri-cube and bi-square functions don't look right to me. Can anyone with more expertise explain this?

1) Tri-cube:
smooth.m line 377 reads:
weight = (1 - (dsort/dmax).^3).^1.5; % tri-cubic weight

Shouldn't that last power be 3 instead of 1.5? Numerous descriptions of the function (including in the Matlab help manual) seem to indicate so.

2) Bi-square:
smooth.m lines 431:439
r1 = abs(r1-median(r1));
mad = median(r1);
if mad > maxabsyXeps
rweight = r1./(6*mad);
id = (rweight<=1);
rweight(~id) = 0;
rweight(id) = (1-rweight(id).*rweight(id));
weight = weight.*rweight;
end

This looks to me like the final "robust weight" factor is simply (1 - (r/(6*MAD))^2), when I believe that entire expression should be squared to match the definition of the bi-square function shown in the Help documentation.

Can anyone clarify this for me? Thanks a lot in advance.
From: John D'Errico on 2 Feb 2010 12:50
"Ryan " <rr2box42(a)gmail.com> wrote in message <hk9n5h$1hf$1(a)fred.mathworks.com>...
> I've got a couple questions about the implementation of loess/lowess in the Matlab smooth function. Specifically, the tri-cube and bi-square functions don't look right to me. Can anyone with more expertise explain this?
>
> 1) Tri-cube:
> smooth.m line 377 reads:
> weight = (1 - (dsort/dmax).^3).^1.5; % tri-cubic weight
>
> Shouldn't that last power be 3 instead of 1.5? Numerous descriptions of the function (including in the Matlab help manual) seem to indicate so.

Often when weights are employed in different schemes,
one must take the sqrt of the weight, depending upon
how the weighting will be implemented.

John
From: Cokelid on 2 Feb 2010 13:22
On Feb 2, 12:25 pm, "Ryan " <rr2bo...(a)gmail.com> wrote:
> I've got a couple questions about the implementation of loess/lowess in the Matlab smooth function.
>
> 1) Tri-cube:
>   smooth.m line 377 reads:  
>     weight = (1 - (dsort/dmax).^3).^1.5; % tri-cubic weight
>
> Shouldn't that last power be 3 instead of 1.5?  Numerous descriptions of the function (including in the Matlab help manual) seem to indicate so.

Hey Ryan,

I don't have the smooth function to check (Curve fitting toolbox?) but
I did implement lowess myself after Cleveland's paper and old Fortran
code. It's possible I missed something but my tri-cubic weights are
calculated as you suggest:
(1 - (abs(x_norm)).^3).^3

I see from the online Matlab documentation that they differentiate
between loess and lowess, does this make any difference?


From: Jon Cherrie on 3 Feb 2010 06:09
Hi Ryan,

As John suggested, the square root of the weights is what is actually required by the implementation. I'm looking at R2009b, and the comment on both the iBisquareWeights and the iTricubeWeights sub-functions inside smooth.m say

% NOTE that this function returns the square root of the weights

On the difference between loess and lowess, this refers to the local polynomial used: lowess for linear, loess for quadratic. This option will not effect the formula used to compute the robust weights.

Cheers,
--Jon

"Ryan " <rr2box42(a)gmail.com> wrote in message <hk9n5h$1hf$1(a)fred.mathworks.com>...
> I've got a couple questions about the implementation of loess/lowess in the Matlab smooth function. Specifically, the tri-cube and bi-square functions don't look right to me. Can anyone with more expertise explain this?
>
> 1) Tri-cube:
> smooth.m line 377 reads:
> weight = (1 - (dsort/dmax).^3).^1.5; % tri-cubic weight
>
> Shouldn't that last power be 3 instead of 1.5? Numerous descriptions of the function (including in the Matlab help manual) seem to indicate so.
>
> 2) Bi-square:
> smooth.m lines 431:439
> r1 = abs(r1-median(r1));
> mad = median(r1);
> if mad > maxabsyXeps
> rweight = r1./(6*mad);
> id = (rweight<=1);
> rweight(~id) = 0;
> rweight(id) = (1-rweight(id).*rweight(id));
> weight = weight.*rweight;
> end
>
> This looks to me like the final "robust weight" factor is simply (1 - (r/(6*MAD))^2), when I believe that entire expression should be squared to match the definition of the bi-square function shown in the Help documentation.
>
> Can anyone clarify this for me? Thanks a lot in advance.
From: Ryan on 3 Feb 2010 09:41
Thanks everyone for the replies. I was in the middle of replying confusedly to Jon because I'm running r2009b too and my smooth.m doesn't have that note nor subfunctions for either weighting algorithm. Then I realized that my smooth.m is about 1,000,000 years old, because a colleague "thoughtfully" packaged one with a program he wrote a few years ago, and I've got that program's directory on the path. So I feel silly about that...

Nevertheless, why do we want to use the square roots of the weights anyway? I've been consulting several references throughout my current project (including Cleveland, 1979 referenced in smooth.m), and I haven't seen such a recommendation so far.

Again, thanks everyone for your time and helpfulness.
 |  Next  |  Last
Pages: 1 2
Prev: Filter and maximal graylevel
Next: [speed up] vectorized find() for sorted data?