AnalysisInterpolationLoess
// This source code is subject to the terms of the Mozilla Public License 2.0 at https://mozilla.org/MPL/2.0/
// © RicardoSantos
// Special apreciation to: @midtownSk8rguy, @kaigouthro
//@version=5

// @description LOESS, local weighted Smoothing function.
library(title='AnalysisInterpolationLoess')

// reference:
// https://commons.apache.org/proper/commons-math/javadocs/api-3.3/org/apache/commons/math3/analysis/interpolation/LoessInterpolator.html

// @function Helper function.
// @param value float, value.
// @returns float.
square (float value) => //{
    math.pow(value, 2.0)
//}

// @function Helper function.
// @param value float, value.
// @returns float.
cube (float value) => //{
    math.pow(value, 3.0)
//}

// @function Helper function for the loess library
// @param value float, value.
// @returns float
tricubic (float value) => //{
    if value > 1.0
        0.0
    else
        cube(1.0 - cube(value))
//}

// @function Helper function for the loess library
// @param value float, value.
// @param exponent float, exponent.
// @returns float
polyweight (float value, float exponent) => //{
    if value > 1.0
        0.0
    else
        math.pow(1.0 - math.pow(value, exponent), exponent)
//{ credit to @midtownSk8rguy
//}}

// @function Helper function for the loess library
// @param value float, value.
// @param exp float, exp.
// @returns float
nweight (float _i1, float _i, float _exp) => //{
    math.exp(-(math.pow(_i1 - _i, _exp)/(math.pow(_exp,_exp))))
// credit to: @kaigouthro
//}


// @function LOESS, local weighted Smoothing function.
// @param sample_x int array, x values.
// @param sample_y float array, y values.
// @param point_span int, local point interval span.
// @returns float array with points y value.
export loess (int[] sample_x, float[] sample_y, int point_span=2) => //{
    int _size_x = array.size(id=sample_x)-1
    int _size_y = array.size(id=sample_y)-1
    // to improve performance there is:
    // no monotonic sequence check, no null checks.
    // Overflow checks:
    switch
        (_size_x < 1 or _size_y < 1) => runtime.error('AnalysisInterpolationLoess -> loess(): samples x and y must have at the least 1 element.')
        (_size_x != _size_y)         => runtime.error('AnalysisInterpolationLoess -> loess(): samples x and y must have the same number of elements.')
        (point_span < 1)             => runtime.error('AnalysisInterpolationLoess -> loess(): point_span must be positive integer number over 0.')
        (2.0 * point_span > _size_x) => runtime.error('AnalysisInterpolationLoess -> loess(): point_span must be smaller than sample size.')
    //
    float[] _yloess = array.new_float(size=1 + _size_x)
    for _i = 0 to _size_x
        int _xi = array.get(id=sample_x, index=_i)
        // define the index interval to pool the values
        int _start = math.max(0, _i - point_span)                               // by: kaigouthro
        int   _end = math.min(_i + point_span, _size_x)                         // by: kaigouthro
        // grab the values for the interval span
        int[]   _span_x = array.slice(id=sample_x, index_from=_start, index_to=_end)
        float[] _span_y = array.slice(id=sample_y, index_from=_start, index_to=_end)
        // calculate the distance, weights then perform local regression
        int _size_s = array.size(id=_span_x) - 1
        float  _sum_w = 0.0
        float  _sum_x = 0.0
        float _sum_x2 = 0.0
        float  _sum_y = 0.0
        float _sum_xy = 0.0
        for _s = 0 to _size_s
            int    _xs = array.get(id=_span_x, index=_s)
            float  _ys = array.get(id=_span_y, index=_s)
            int    _ds = math.abs(_xs-_xi)
            float  _ws = tricubic(_ds / point_span)
            float _wxs = _ws * _xs
            _sum_w  += _ws
            _sum_x  += _wxs
            _sum_x2 += _wxs * _xs
            _sum_y  += _ws  * _ys
            _sum_xy += _wxs * _ys
        float  _mean_x = _sum_x  / _sum_w
        float  _mean_y = _sum_y  / _sum_w
        float _mean_xy = _sum_xy / _sum_w
        float _mean_x2 = _sum_x2 / _sum_w
        float   _mx2xx = _mean_x2 - square(_mean_x)
        float    _beta = 0.0
        if math.sqrt(math.abs(_mx2xx)) < 1.0e-15
            _beta := (_mean_xy - _mean_x * _mean_y) / _mx2xx
        float _alpha = _mean_y - _beta   * _mean_x
        //update current y loess value at _i index
        array.set(id=_yloess, index=_i, value=_beta * _xi + _alpha)
    _yloess
//{ usage:
// length = input.int(50)
// npoints= input.int(7)
// f_l1(_x1, _y1, _x2, _y2) => line.new(_x1, _y1, _x2, _y2, color=#DDDD00, style=line.style_solid, width=1)

// import RicardoSantos/ArrayExtension/2 as ae
// a = bar_index < length ? array.from(0) : ae.to_int(ae.sequence_from_series(bar_index, length, 0, true), 'floor')
// b = bar_index < length ? array.from(0.0) : ae.sequence_from_series(close, length, 0, true)

// if barstate.islastconfirmedhistory
//     c = loess(a, b, npoints)
//     label.new(bar_index, 0.0, str.tostring(c, '#.###'))


// 	for _i = 1 to array.size(a)-1
//         _x1 =  array.get(a, _i-1)
//         _x2 =  array.get(a, _i  )
//         _by1 = array.get(c, _i-1)
//         _by2 = array.get(c, _i)
//         f_l1(_x1, _by1, _x2, _by2)

// plot(close)
//}}

// @function aLOESS, adaptive local weighted Smoothing function.
// @param sample_x int array, x values.
// @param sample_y float array, y values.
// @param point_span int, local point interval span.
// @param span_exp float, local span exp factor, default=na, will run regular LOESS algorithm, adaptive otherwise.
// @returns float array with points y value.
export aloess (int[] sample_x, float[] sample_y, int point_span=2, float span_exponent=3.0) => //{
    int _size_x = array.size(id=sample_x)-1
    int _size_y = array.size(id=sample_y)-1
    // to improve performance there is:
    // no monotonic sequence check, no null checks.
    // Overflow checks:
    switch
        (_size_x < 1 or _size_y < 1) => runtime.error('AnalysisInterpolationLoess -> loess(): samples x and y must have at the least 1 element.')
        (_size_x != _size_y)         => runtime.error('AnalysisInterpolationLoess -> loess(): samples x and y must have the same number of elements.')
        (point_span < 1)             => runtime.error('AnalysisInterpolationLoess -> loess(): point_span must be positive integer number over 0.')
        (2.0 * point_span > _size_x) => runtime.error('AnalysisInterpolationLoess -> loess(): point_span must be smaller than sample size.')
    //
    float[] _yloess = array.new_float(size=1 + _size_x)
    for _i = 0 to _size_x
        int _xi = array.get(id=sample_x, index=_i)
        // define the index interval to pool the values
        int _start = math.max(0, _i - point_span)                               // by: kaigouthro
        int   _end = math.min(_i + point_span, _size_x)                         // by: kaigouthro
        // grab the values for the interval span
        int[]   _span_x = array.slice(id=sample_x, index_from=_start, index_to=_end)
        float[] _span_y = array.slice(id=sample_y, index_from=_start, index_to=_end)
        // calculate the distance, weights then perform local regression
        int   _size_s = array.size(id=_span_x)-1
        float  _sum_w = 0.0
        float  _sum_x = 0.0
        float _sum_x2 = 0.0
        float  _sum_y = 0.0
        float _sum_xy = 0.0
        for _s = 0 to _size_s
            int    _xs = array.get(id=_span_x, index=_s)
            float  _ys = array.get(id=_span_y, index=_s)
            int    _ds = math.abs(_xs - _xi)
            float  _ws = polyweight(_ds / point_span, span_exponent)
            float _wxs = _ws * _xs
            _sum_w  += _ws
            _sum_x  += _wxs
            _sum_x2 += _wxs * _xs
            _sum_y  += _ws  * _ys
            _sum_xy += _wxs * _ys
        float  _mean_x = _sum_x  / _sum_w
        float  _mean_y = _sum_y  / _sum_w
        float _mean_xy = _sum_xy / _sum_w
        float _mean_x2 = _sum_x2 / _sum_w
        float   _mx2xx = _mean_x2 - square(_mean_x)
        float    _beta = 0.0
        if math.sqrt(math.abs(_mx2xx)) < 1.0e-15
            _beta := (_mean_xy - _mean_x * _mean_y) / _mx2xx
        float _alpha = _mean_y - _beta   * _mean_x
        //update current y loess value at _i index
        array.set(id=_yloess, index=_i, value=_beta * _xi + _alpha)
    _yloess
//{ usage:
adaptive = input.bool(false)
length = input.int(50)
npoints = input.int(7)
exp = input.float(3.0)
f_l1(_x1, _y1, _x2, _y2) => line.new(_x1, _y1, _x2, _y2, color=#DDDD00, style=line.style_solid, width=1)

import RicardoSantos/ArrayGenerate/1 as ag
import RicardoSantos/ArrayExtension/1 as ae

a = bar_index <= length ? array.from(0) : ae.to_int(ag.sequence_from_series(bar_index, length, 0, true), 'floor')
b = bar_index <= length ? array.from(0.0) : ag.sequence_from_series(close, length, 0, true)
if barstate.islastconfirmedhistory
    c = adaptive ? aloess(a, b, npoints, exp) : loess(a, b, npoints)
    // label.new(bar_index, 0.0, str.tostring(c, '#.###'))


	for _i = 1 to array.size(a)-1
        _x1 =  array.get(a, _i-1)
        _x2 =  array.get(a, _i  )
        _by1 = array.get(c, _i-1)
        _by2 = array.get(c, _i)
        f_l1(_x1, _by1, _x2, _by2)

plot(close)
//}}

// @function LOESS moving average.
// @param source float, data source.
// @param length int, data source sampling length.
// @param point_span int, local point interval span.
// @param span_exp float, local span exp factor, default=na, will run regular LOESS algorithm, adaptive otherwise.
// @retuns float.
// Note: this only pools the last value of the series so it will not produce the smoothed results as in the original.
ma (float source, int length, int point_span=3, float span_exp=na) => //{
    var int[]   _x = array.new_int(length, bar_index)
    var float[] _y = array.new_float(length, source)
    array.pop(_x), array.unshift(_x, bar_index)
    array.pop(_y), array.unshift(_y, source)
    if na(span_exp)
        array.get(loess(_x, _y, point_span), 0)
    else
        array.get(aloess(_x, _y, point_span, span_exp), 0)
//{ usage:
plot(ma(close, length, npoints, adaptive?exp:na), color=color.red, linewidth=3)
//}}

// @function LOESS moving average by kaigouthro.
// @param source float, data source.
// @param len int, data source sampling length.
// @param point_spaninput int, local point interval span.
// @param exp float, local span exp factor, default=0.0, will run regular LOESS algorithm, adaptive otherwise.
// @retuns float.
loest (float src, int len=18, int point_spaninput=15, float exp=0.0) => //{
    float[] _out_array =  array.new_float()
    _length    =  len - (point_spaninput)
    point_span =  math.max(2, math.min(math.abs(point_spaninput), math.floor(len/3)))
    tlen   = math.min(_length, len-point_span)
    for _i = 0 to tlen by 1
        _sum_w  = 0., _sum_x = 0., _sum_x2 = 0., _sum_y = 0., _sum_xy = 0., _beta = 0.
        int _end    = math.min(_i + point_span, len - 1)
        int _start  = math.max(0, _end - point_span)
        int _xi     = tlen - _i
        for _s = 0 to point_span -1
            int    _xs = tlen - (_i + _s)
            _ys        = src[_i +_s]
            _ds        = math.abs( _xs  - _xi) *  (point_span / _length) * math.pow(math.phi, math.phi)
            _ws        = nweight(_ds, point_span, exp)
            _wxs       = _ws * _xs * (exp * 1.5-((1 + _s)/point_span))
            _sum_w    += _ws
            _sum_x    += _wxs              , _sum_x2 += _wxs * _xs    
            _sum_y    += _ws      * _ys    , _sum_xy += _wxs * _ys    
        float _mean_x  = _sum_x   / _sum_w 
        float _mean_y  = _sum_y   / _sum_w 
        float _mean_xy = _sum_xy  / _sum_w 
        float _mean_x2 = _sum_x2  / _sum_w 
        float _mx2xx   = _mean_x2 - _mean_x * _mean_x 
        if math.sqrt(math.abs(_mx2xx)) < 1.0e-15
            _beta := (_mean_xy - _mean_x * _mean_y) / _mx2xx
        float _alpha = _mean_y - _beta   * _mean_x
        array.push(_out_array, math.round_to_mintick(_beta * _xi + _alpha))
    array.avg(_out_array)
// usage:
oc  = loest(close, length, npoints, adaptive?exp:0.0)
plot(oc, 'loesst', oc>oc[1] ? color.green:color.orange,2)
// credit to: @kaigouthro.
//}