Python smoothers_lowess.lowess函数代码示例

def plotting(file, ext):
    ts = pd.read_csv(file, header=[0, 1], index_col=[0])
    ts.index = [dt.datetime.strptime(x, "%H:%M:%S").time() for x in ts.index]
    ts.columns.set_levels(
        [dt.datetime.strptime(x, "%Y-%m-%d").date() for x in ts.columns.levels[0].values], 0, inplace=True
    )

    # gets used in the smoothing
    xtime = [int(x.hour) + int(x.minute) / 60 for x in ts.index]
    timeind = pd.date_range("00:00", "23:59", freq="min").to_pydatetime()  ###Very important!

    # for each month in the data, construct a df of just that month
    for mo in set([x.month for x in ts.columns.levels[0]]):
        monthnum = mo
        motxt = time.strftime("%B", time.strptime(str(monthnum), "%m"))
        mocheck = [x.month == monthnum for x in ts.columns.levels[0]]
        moflag = []
        for ans in mocheck:
            moflag.append(ans)
            moflag.append(ans)

        onemo = ts.loc[:, moflag]

        sleepy = onemo.xs("Sleep", level=1, axis=1).sum(axis=1) / (len(onemo.columns) / 2)
        eaty = onemo.xs("Eat", level=1, axis=1).sum(axis=1) / (len(onemo.columns) / 2)

        # begin plotting
        fig = plt.figure(figsize=(18, 6))
        ax = fig.add_subplot(111)

        ####Plot Sleep
        filtereds = lowess(sleepy, xtime, is_sorted=True, frac=0.025, it=0)
        ax.plot(timeind, filtereds[:, 1], "b", linewidth=2, label="Sleeping")
        ax.fill_between(timeind, 0, filtereds[:, 1], alpha=0.3, facecolor="b")
        # ax.plot(ts.index,sleepy,'b',linewidth=2,label='Sleeping')#raw data, not smoothed
        # ax.fill_between(ts.index, 0, sleepy,alpha=0.3,facecolor='b')

        ####Plot Eat
        filterede = lowess(eaty, xtime, is_sorted=True, frac=0.025, it=0)
        ax.plot(timeind, filterede[:, 1], "orange", linewidth=2, label="Eating")
        ax.fill_between(timeind, 0, filterede[:, 1], alpha=0.3, facecolor="orange")
        # ax.plot(ts.index,eaty,'orange',linewidth=2,label='Eating')
        # ax.fill_between(ts.index, 0, eaty,alpha=0.3,facecolor='orange')

        ####Axis formatting
        xax = ax.get_xaxis()
        xax.set_major_locator(mdates.HourLocator(byhour=range(0, 24, 2)))
        xax.set_major_formatter(mdates.DateFormatter("%H:%M"))
        ax.set_title("Activity Fraction at a Given Time of Day", fontsize="xx-large")
        ax.text("16:00", 0.9, motxt, fontsize="xx-large", color="k", fontweight="bold")
        ax.legend(fontsize="x-large")
        ax.set_ylim(0, 1.1)
        fig.autofmt_xdate()
        filename = "b2_TimeSeries/Activity_" + str(monthnum) + "." + ext
        fig.savefig(filename)
    return

def loess(x,y,frac=0.2,it=None,scatter=True):
    from statsmodels.nonparametric.smoothers_lowess import lowess 
    y = np.array(y)
    x = np.array(x)
    y = y[x.argsort()] # Sort y according to order of x.  
    x.sort() # Sort x in place.  
    if it is not None: # Helps if you are getting NaN's in the output.  
        d = lowess(y,x,frac=frac,it=it)
    else:
        d = lowess(y,x,frac=frac)
    return d

    def test_frac(self):
        rfile = os.path.join(rpath, "test_lowess_frac.csv")
        test_data = np.genfromtxt(open(rfile, "rb"), delimiter=",", names=True)

        expected_lowess_23 = np.array([test_data["x"], test_data["out_2_3"]]).T
        expected_lowess_15 = np.array([test_data["x"], test_data["out_1_5"]]).T

        actual_lowess_23 = lowess(test_data["y"], test_data["x"], frac=2.0 / 3)
        actual_lowess_15 = lowess(test_data["y"], test_data["x"], frac=1.0 / 5)

        assert_almost_equal(expected_lowess_23, actual_lowess_23, decimal=testdec - 1)
        assert_almost_equal(expected_lowess_15, actual_lowess_15, decimal=testdec)

    def test_iter(self):
        rfile = os.path.join(rpath, "test_lowess_iter.csv")
        test_data = np.genfromtxt(open(rfile, "rb"), delimiter=",", names=True)

        expected_lowess_no_iter = np.array([test_data["x"], test_data["out_0"]]).T
        expected_lowess_3_iter = np.array([test_data["x"], test_data["out_3"]]).T

        actual_lowess_no_iter = lowess(test_data["y"], test_data["x"], it=0)
        actual_lowess_3_iter = lowess(test_data["y"], test_data["x"], it=3)

        assert_almost_equal(expected_lowess_no_iter, actual_lowess_no_iter, decimal=testdec)
        assert_almost_equal(expected_lowess_3_iter, actual_lowess_3_iter, decimal=testdec)

    def test_frac(self):
        rfile = os.path.join(rpath, 'test_lowess_frac.csv')
        test_data = np.genfromtxt(open(rfile, 'rb'),
                                  delimiter = ',', names = True)

        expected_lowess_23 = np.array([test_data['x'], test_data['out_2_3']]).T
        expected_lowess_15 = np.array([test_data['x'], test_data['out_1_5']]).T

        actual_lowess_23 = lowess(test_data['y'], test_data['x'] ,frac = 2./3)
        actual_lowess_15 = lowess(test_data['y'], test_data['x'] ,frac = 1./5)

        assert_almost_equal(expected_lowess_23, actual_lowess_23, decimal = testdec-1)
        assert_almost_equal(expected_lowess_15, actual_lowess_15, decimal = testdec)

    def test_iter(self):
        rfile = os.path.join(rpath, 'test_lowess_iter.csv')
        test_data = np.genfromtxt(open(rfile, 'rb'),
                                  delimiter = ',', names = True)

        expected_lowess_no_iter = np.array([test_data['x'], test_data['out_0']]).T
        expected_lowess_3_iter = np.array([test_data['x'], test_data['out_3']]).T

        actual_lowess_no_iter = lowess(test_data['y'], test_data['x'], it = 0)
        actual_lowess_3_iter = lowess(test_data['y'], test_data['x'], it = 3)

        assert_almost_equal(expected_lowess_no_iter, actual_lowess_no_iter, decimal = testdec)
        assert_almost_equal(expected_lowess_3_iter, actual_lowess_3_iter, decimal = testdec)

    def test_delta(self):
        rfile = os.path.join(rpath, "test_lowess_delta.csv")
        test_data = np.genfromtxt(open(rfile, "rb"), delimiter=",", names=True)

        expected_lowess_del0 = np.array([test_data["x"], test_data["out_0"]]).T
        expected_lowess_delRdef = np.array([test_data["x"], test_data["out_Rdef"]]).T
        expected_lowess_del1 = np.array([test_data["x"], test_data["out_1"]]).T

        actual_lowess_del0 = lowess(test_data["y"], test_data["x"], frac=0.1)
        actual_lowess_delRdef = lowess(test_data["y"], test_data["x"], frac=0.1, delta=0.01 * np.ptp(test_data["x"]))
        actual_lowess_del1 = lowess(test_data["y"], test_data["x"], frac=0.1, delta=1.0 + 1e-10)

        assert_almost_equal(expected_lowess_del0, actual_lowess_del0, decimal=testdec)
        assert_almost_equal(expected_lowess_delRdef, actual_lowess_delRdef, decimal=testdec)
        assert_almost_equal(expected_lowess_del1, actual_lowess_del1, decimal=10)  # testdec)

    def test_simple(self):
        x = np.arange(20, dtype='float32')
        #standard normal noise
        noise = np.array([-0.76741118, -0.30754369,
                            0.39950921, -0.46352422, -1.67081778,
                            0.6595567 ,  0.66367639, -2.04388585,
                            0.8123281 ,  1.45977518,
                            1.21428038,  1.29296866,  0.78028477,
                            -0.2402853 , -0.21721302,
                            0.24549405,  0.25987014, -0.90709034,
                            -1.45688216, -0.31780505])
        y = x + noise

        # R output
        out = [-0.6260344553, 0.565071712, 1.759627189,
                2.9579633258, 4.1560636154, 5.3473396937,
                6.522298218, 7.708159388, 8.8759055519,
                9.9409758603, 10.8981138458, 11.7851424728,
                12.6188717297, 13.4098497374, 14.1516996585,
                14.9180658147, 15.6956600199, 16.4783034134,
                17.2617441531, 18.0459201716]

        expected_lowess = np.array([x, out]).T

        actual_lowess = lowess(y,x)

        assert_almost_equal(expected_lowess, actual_lowess)

def add_lowess(ax, lines_idx=0, frac=.2, **lowess_kwargs):
    """
    Add Lowess line to a plot.

    Parameters
    ----------
    ax : matplotlib Axes instance
        The Axes to which to add the plot
    lines_idx : int
        This is the line on the existing plot to which you want to add
        a smoothed lowess line.
    frac : float
        The fraction of the points to use when doing the lowess fit.
    lowess_kwargs
        Additional keyword arguments are passes to lowess.

    Returns
    -------
    fig : matplotlib Figure instance
        The figure that holds the instance.
    """
    y0 = ax.get_lines()[lines_idx]._y
    x0 = ax.get_lines()[lines_idx]._x
    lres = lowess(y0, x0, frac=frac, **lowess_kwargs)
    ax.plot(lres[:, 0], lres[:, 1], 'r', lw=1.5)
    return ax.figure

def regression_plot(Z,X,band_names=None,visible_only=True,figsize=(12,7)):
    """
    Produce a figure with a plot for each image band that displays the
    relationship between depth and radiance and gives a visual representation
    of the regression carried out in the `slopes` and `regressions` methods.

    Notes
    -----
    This method doesn't come directly from Lyzenga 1978 but the author of this
    code found it helpful.

    Parameters
    ----------
    Z : np.ma.MaskedArray
        Array of depth values repeated for each band so that Z.shape==X.shape.
        The mask needs to be the same too so that Z.mask==X.mask for all the
        bands.
    X : np.ma.MaskedArray
        The array of log transformed radiance values from equation B1 of
        Lyzenga 1978.

    Returns
    -------
    figure
        A matplotlib figure.
    """
    if band_names is None:
        band_names = ['Band'+str(i+1) for i in range(X.shape[-1])]
    nbands = X.shape[-1]
    if np.atleast_3d(Z).shape[-1] == 1:
        Z = np.repeat(np.atleast_3d(Z), nbands, 2)
    if visible_only:
        fig, axs = plt.subplots( 2, 3, figsize=figsize)
    else:
        fig, axs = plt.subplots( 2, 4, figsize=figsize )
    regs = regressions(Z,X)
    for i, ax in enumerate(axs.flatten()):
        if i > nbands-1:
            continue
        slp, incpt, rval = regs[:,i]
        # print X.shape, Z.shape
        x, y = equalize_array_masks(Z[...,i], X[...,i])
        if x.count() < 2:
            continue
        x, y = x.compressed(), y.compressed()
        # print "i = {}, x.shape = {}, y.shape = {}".format(i, x.shape, y.shape)
        ax.scatter( x, y, alpha=0.1, edgecolor='none', c='gold' )
        smth = lowess(y,x,frac=0.2)
        # ax.plot(smth.T[0],smth.T[1],c='black',alpha=0.5)
        ax.plot(smth.T[0],smth.T[1],c='black',alpha=0.5,linestyle='--')
        reglabel = "m=%.2f, r=%.2f" % (slp,rval)
        f = lambda x: incpt + slp * x
        ax.plot( x, f(x), c='brown', label=reglabel, alpha=1.0 )
        ax.set_title( band_names[i] )
        ax.set_xlabel( r'Depth (m)' )
        ax.set_ylabel( r'$X_i$' )
        ax.legend(fancybox=True, framealpha=0.5)
    plt.tight_layout()
    return fig

    def test_delta(self):
        rfile = os.path.join(rpath, 'test_lowess_delta.csv')
        test_data = np.genfromtxt(open(rfile, 'rb'),
                                  delimiter = ',', names = True)

        expected_lowess_del0 = np.array([test_data['x'], test_data['out_0']]).T
        expected_lowess_delRdef = np.array([test_data['x'], test_data['out_Rdef']]).T
        expected_lowess_del1 = np.array([test_data['x'], test_data['out_1']]).T

        actual_lowess_del0    = lowess(test_data['y'], test_data['x'], frac=0.1)
        actual_lowess_delRdef = lowess(test_data['y'], test_data['x'], frac=0.1,
                       delta = 0.01 * np.ptp(test_data['x']))
        actual_lowess_del1    = lowess(test_data['y'], test_data['x'], frac = 0.1, delta = 1.0 + 1e-10)

        assert_almost_equal(expected_lowess_del0, actual_lowess_del0, decimal = testdec)
        assert_almost_equal(expected_lowess_delRdef, actual_lowess_delRdef, decimal = testdec)
        assert_almost_equal(expected_lowess_del1, actual_lowess_del1, decimal = 10) #testdec)

    def test_simple(self):
        rfile = os.path.join(rpath, "test_lowess_simple.csv")
        test_data = np.genfromtxt(open(rfile, "rb"), delimiter=",", names=True)
        expected_lowess = np.array([test_data["x"], test_data["out"]]).T

        actual_lowess = lowess(test_data["y"], test_data["x"])

        assert_almost_equal(expected_lowess, actual_lowess, decimal=testdec)

    def test_simple(self):
        rfile = os.path.join(rpath, 'test_lowess_simple.csv')
        test_data = np.genfromtxt(open(rfile, 'rb'),
                                  delimiter = ',', names = True)
        expected_lowess = np.array([test_data['x'], test_data['out']]).T

        actual_lowess = lowess(test_data['y'], test_data['x'])

        assert_almost_equal(expected_lowess, actual_lowess, decimal = testdec)

 def generate(name, fname, x='x', y='y', out='out', kwargs=None, decimal=7):
     kwargs = {} if kwargs is None else kwargs
     data = np.genfromtxt(os.path.join(rpath, fname), delimiter=',', names=True)
     assert_almost_equal.description = name
     if callable(kwargs):
         kwargs = kwargs(data)
     result = lowess(data[y], data[x], **kwargs)
     expect = np.array([data[x], data[out]]).T
     assert_almost_equal(result, expect, decimal)

 def generate(name, fname,
              x='x', y='y', out='out', kwargs={}, decimal=7):
     data = np.genfromtxt(
         os.path.join(rpath, fname), delimiter=',', names=True)
     assert_equal_at_testdec = partial(
         assert_almost_equal, decimal=decimal)
     assert_equal_at_testdec.description = name
     if callable(kwargs):
         kwargs = kwargs(data)
     result = lowess(data[y], data[x], **kwargs)
     expect = np.array([data[x], data[out]]).T
     return assert_equal_at_testdec, result, expect