Bookmark and Share

Notice: On March 31, it was announced that Statalist is moving from an email list to a forum. The old list will shut down at the end of May, and its replacement, statalist.org is already up and running.


[Date Prev][Date Next][Thread Prev][Thread Next][Date Index][Thread Index]

st: trying to understand -lpoly- with degree zero


From   Dan Weitzenfeld <dan.weitzenfeld@emsense.com>
To   "statalist@hsphsun2.harvard.edu" <statalist@hsphsun2.harvard.edu>
Subject   st: trying to understand -lpoly- with degree zero
Date   Mon, 12 Apr 2010 18:19:56 -0500

Hi All,
I'm trying to reverse-engineer the lpoly command, so that I can code it in a different language.
Specifically, I'm trying to mimic the command with the default options of degree(0) and kernel(epanechnikov).
Within -help lpoly-, it says:


    degree(#) specifies the degree of the polynomial to be used in the smoothing.  The default is degree(0), meaning local-mean
        smoothing.

I understand this to mean that at each point in the smoothing grid, a kernel-weighted mean is taken.  For my purposes, it makes sense to use every point in the dataset as a point in the grid.  But when I tried to code that - see below - my resulting curve is less smooth than the output of -lpoly- with degree 0, same kernel, same half-width.

According to wikipedia, the epanechnikov kernel is:
K(u) = (3/4)(1-u^2) 
Am I scaling the kernel incorrectly?

Thanks in advance, 
Dan

*----------------------------------------------------

webuse motorcycle, clear
*my dataset is a constant 60 hz; changing this dataset to mimic
replace time = _n*(1/60)
local  halfwidth = .04
local  halfwidth_squared = `halfwidth'*`halfwidth'
qui gen fake_lpoly = .
qui count
local vals = r(N)
forvalues i = 1(1)`vals' {
 local time = time[`i']
 capture drop distance* kernel product
 qui gen distance_from_point = time-`time'
 qui gen distance_sq =  distance_from_point*distance_from_point
*generating epanechnikov kernel, scaled to local halfwidth
 qui gen kernel = (3/4)*(1- (distance_sq/`halfwidth_squared'))
 qui gen product = kernel*accel if kernel>=0
 qui sum product
 local numerator = r(sum)
 qui sum kernel if kernel>=0
 local denominator = r(sum)
 qui replace fake_lpoly = `numerator'/`denominator' in `i'
}

*compare with -LPOLY-
sort time
capture drop x real_lpoly
lpoly accel time, bw(`halfwidth') n(`vals') gen(x real_lpoly) nogra
twoway (line  real_lpoly  fake_lpoly time) (scatter  accel time)

*--------------------------------------------------

This email message is confidential and is intended only for the named recipient(s) above, and may contain information that is privileged, attorney work product or exempt from disclosure under applicable law.  The contents of this email are for discussion purposes only and shall not create a legally binding commitment or agreement of EmSense Corporation, which shall require EmSense Corporation management approval and a formal written agreement signed by both parties.  If you have received this message in error, or are not the named recipient(s), please immediately notify the sender and delete this message from your computer.   

*
*   For searches and help try:
*   http://www.stata.com/help.cgi?search
*   http://www.stata.com/support/statalist/faq
*   http://www.ats.ucla.edu/stat/stata/


© Copyright 1996–2014 StataCorp LP   |   Terms of use   |   Privacy   |   Contact us   |   Site index