# Kernel density regression in Julia

Here is a naive implementation of a kernel density regression algorithm in Julia.

The code works, but there is quite a bit of room for improvement though. First, there is not an optimal bandwidth selection implemented and second; it's really slow. Nevertheless, it took me a while to implement and I am really looking forward to any suggestions on how to improve the precision and performance of the code.

function kernelRegression(x, y, h)
#length of the actual data
n = length( x )

#if length of the data is smaller then 100,
#use 100 for the length of the sequence to make sure the curve
#can be plotted smoothly
l = min( n, 100 )
seq = linspace(minimum(x), maximum(x), l)

#bandwidth for the gaussioan kernel approx to:
h = 1.06 * std( x ) * n ^(-.2)

#store the result in
est = zeros(Real, l)

#gaussian kernel
function kernel( u )
exp ( -.5 * u^2 )
end

#compute the estimate over the entire sequence of values
#basic algorithm:
#for every x in sequence from min ( data ) to max ( data )
#   for every observation in data
#       calulate the kernel density function at this particular point
#       which is the sum of all the kernels on the x axis devided by n * binwidth
#       multiply the kernel at this point with the respective y value
#       and devide by the sum of kernels
for (index, sequence) in enumerate( seq )

#reset the weights
weight = 0.0

for ( i, ob ) in enumerate( x )

temp = 0.0
#inner loop: calculate the mean
#of the kernels at X_i = x
for obs in  x
u = ( sequence - obs ) / h
temp += kernel ( u )
end
temp /= ( n * h )

#use temp to calculate the estimate of y at X_i = x
v = ( sequence - ob ) / h
weight += ( kernel ( v ) / temp ) * y[i]
end

est[index] = weight / ( n * h )
end
est
#scatter(x, y)
#plot(seq, est, color="red")
end

(() ->
begin
x = rand(Uniform(-6, 6), 1000)
y = cos(x) + rand(Normal(0, 0.2), 1000)
kernelRegression(x, y)
scatter(x,y)
plot(linspace(-6,6,1000), y_est)
end)()