本文介绍了python中的完全单调插值的处理方法,对大家解决问题具有一定的参考价值,需要的朋友们下面随着小编来一起学习吧!

问题描述

我有一些数据,例如:

我想为它拟合一条可微的单调曲线.我尝试了

这不是单调的.

如何为这样的数据拟合单调曲线?

这是另一个类似图形的 y 值样本集:

[0.1109157119023644,0.20187393816931934,0.14466318670239758,0.16535159414166822,0.05452708697483864,0.2153046237959556,0.2200300476272603,0.21012762463269324,0.15947100322395022,0.2819691842129948,0.15567770052985092,0.24850595803020692,0.1329341593280457,0.15595107081606913,0.3232021121832229,0.23707961921686588,0.2415887076540357,0.32363506549779797,0.3584089204036798,0.29232772580068433,0.22145994836140775,0.22797587985241133,0.2717787840603025,0.3245255944762287,0.29301098282789195,0.32417076823344143,0.3450906550996232,0.34272097408024904,0.3868714875012437,0.41876692320045755,0.3544198724867363,0.33073960954801895,0.3921033666371904,0.33349050060172974,0.3608862044547096,0.37375822841635425,0.5396399750708429,0.4209201143798284,0.42004773793166883,0.5217725632679073,0.5911731474218788,0.43389609315065386,0.4287288396176006,0.43007525393257007,0.5687062142675405,0.6030811498722173,0.5292225577714743, 0.47710974351051355, 0.6182720730381119,0.6241033581931327,0.6236788197617511,0.6643161356364049,0.5577616524049582,0.6888440258481371,0.6867893120660341,0.6685257606057502,0.599481675493677,0.7309075091448749,0.7644365338580481,0.6176797601816733,0.6751467827192018,0.6452178017908761,0.6684778262246701,0.7003380077556168,0.667035916425416,0.8434451759113093,0.8419343615815968,0.8657695361433773,0.7392487161484605,0.8773282098364621,0.8265679895117846,0.7246599961191632,0.7251899061730714,0.9271640780410231,0.9180581424305536,0.8099033021701689,0.8268585329594615,0.8519967080830176,0.8711231413093845,0.8689802343798663,0.8299523829217353,1.0057741699770046,0.8538130788729608,0.9662784297225102,1.023419780920539,0.913146849759822,0.9900885996579213,0.8740638988529978,0.8900285618419457,0.9065474574434158,1.0749522597307315,1.0319120938258166,1.0051369663172995,0.9893558841613622,1.051384986916457,1.0327996870915341,1.0945543972861898,0.9716604944496021,1.1490370559566179,1.1379231481207432,1.6836433783615088,1.8162068766097395,2.072155286917785,2.0395966998366,2.191064589600466,2.1581974932543617,2.163403843819597,2.133441151300847,2.1726053994136922,2.1157865673629526,2.2249636455682866,2.2313062166802147,2.1731708496472764,2.315203950110816,2.1601242661726827,2.174940281421225,2.2653635413275945,2.337227057574145,2.3645767548381618,2.3084919291392527,2.314014515926446,2.25166717296155,2.2621157708115778,2.2644578546265586,2.313504860292943,2.398969190357051,2.309443951779675,2.278946047410807,2.4080802287121146,2.353652872018618,2.35527529074088,2.4233001060410784,2.428767198055608,2.35677123091093,2.497135132404064,2.3978099128437282,2.3970802609341972,2.4967434818740024,2.511209192435555,2.541001050440798,2.5760248002036525,2.5960512284192245,2.4778408861721037,2.5757724103530046,2.631148267999664,2.538327346218921,2.4878734713248507,2.6133797275761066,2.6282561527857395、2.6150327104952447,3.102757164382848,3.3318503012160905,3.3907776288198193,3.6065313558941936,3.601180295875859,3.560491539319038,3.650095006265445,3.574812155815713,3.686227315374108,3.6338261415040867,3.5661194785086288,3.5747332336054645,3.560674343726918,3.5678550481603635,3.5342848534390967,3.4929538312485913,3.564544653619436,3.6861775399566126,3.6390300636595216,3.6656336332413666,3.5731185631923945,3.5965520044069854,3.537434489989021,3.5590937423870144,3.5331656424410083,3.640652819618705,3.5971240740252126,3.641793843012055,3.6064014089254295,3.530378938786505,3.613631139461306,3.519542268056021,3.5416251524576,3.524789618934195,3.5519951806099512,3.6435695455293975,3.6825670484650863,3.5993379768209217,3.628367553897596,3.633290480934276,3.5772841681579535,3.602326323397947,3.518180278272883,3.531054006706696,3.5566645495066167,3.5410992153240985,3.630762839301216,3.5924649123201053,3.646230633817883,3.568290612034935、3.638356129262967、3.566083243271712,3.6064978645771797,3.4942864293427633,3.595438454812999,3.681726879126678,3.6501308156903463,3.5490717955938593,3.598535359345363,3.6328331698421654,3.595159538698094,3.556715819008055,3.6292942886764554,3.6362895697392856,3.5965220100874093,3.6103542985016266,3.5715010140382493,3.658769915445062,3.5939686395400416,3.4974461928859917,3.5232691556732267,3.6145687814416614,3.5682054018341005,3.648937250575395,3.4912089018613384,3.522426560340423,3.6757968409374637,3.651348691084845,3.5395070091675973,3.5306275536360383,3.6153498246329883,3.599762785949876,3.5351931286962333,3.6488316987683054,3.5198301490992963,3.5696570079786687,3.561553836008927,3.5659475947331423,3.553147100256108,3.5475591872743664,3.6097226797553317,3.6849600324757934,3.5264731043844413,3.506658609738451,3.5535775980874114,3.5487291053913554,3.570651383823912,3.552993371839188,3.5054297764661846,3.5723024888238792]
解决方案

这里有一个基本上是 5 行 python 的单调曲线拟合器,带有 numpy和来自 scipy.signal 的低通滤波器:

#!/usr/bin/env python"""https://stackoverflow.com/questions/56551114/fully-monotone-interpolation-in-python"""# 也可以看看# https://en.wikipedia.org/wiki/Monotone-spline 又名 I-spline# https://scikit-learn.org/stable/modules/isotonic.html# 丹尼斯 2020 年 3 月 2 日从 __future__ 导入师,print_function将 numpy 导入为 np从 scipy 导入信号作为 sig从 matplotlib 导入 pyplot 作为 plt将 seaborn 作为 sns 导入def butter_filtfilt( x, Wn=0.5, axis=0 ):""" 黄油 ( 2, Wn ), filtfilt轴 0 每列,-1 每行"""b, a = sig.butter( N=2, Wn=Wn )return sig.filtfilt( b, a, x, axis=axis, method="gust" ) # 两次,向前向后定义整数( x ):返回 x.round().astype(int)def minavmax( x ):返回最小 av 最大 %.3g %.3g %.3g" % (x.min(), x.mean(), x.max() )def pvec( x ):n = len(x)//25 * 25返回 "%s \n%s \n" % (minavmax( x ),ints( x[ - n : ]) .reshape( -1, 25 ))#………………………………………………………………………………………………......................................def monofit( y, Wn=0.1 ):"""单调递增曲线拟合"""y = np.asarray(y).squeeze()打印(\n{ monofit:y %d %s Wn %.3g" % (len(y), minavmax( y ), Wn ))ygrad = np.gradient( y )打印(grad y:",pvec(ygrad))# 低通滤波器  -gradsmooth = butter_filtfilt( ygrad, Wn=Wn )打印(梯度平滑:",pvec(梯度平滑))ge0 = np.fmax(梯度平滑,0)ymono = np.cumsum( ge0 ) # 整合,对前几个敏感ymono += (y - ymono).mean()错误 = y - ymono打印(y - ymono:",pvec(错误))errstr = "平均值 |y - monofit|: %.2g" % np.abs( err ).mean()打印( errstr )打印(}\n")返回 ymono, err, errstr#………………………………………………………………………………………………......................................如果 __name__ == "__main__":导入系统np.set_printoptions( threshold=20, edgeitems=15, linewidth=120,formatter = dict( float = lambda x: "%.2g" % x )) # 浮点数组 %.2g打印(80 *=")thispy = sys.argv[0]infile = sys.argv[1] 如果 len(sys.argv) >1 \否则so-mono.txt"Wn = 0.1 # ?params = "%s %s Wn %g" % (thispy, infile, Wn)打印(参数)y = np.loadtxt( infile ) * 100打印(y:",y)ymono, err, errstr = monofit( y, Wn=Wn )如果 1:sns.set_style("whitegrid")图, ax = plt.subplots( figsize=[10, 5] )plt.subplots_adjust( left=.05, right=.99, bottom=.05, top=.90 )fig.suptitle(简单的单调曲线拟合:np.gradient | 低通滤波器 | 剪辑  y )[0]plt.vlines( xax[j], y[j], ymono[j], color="blue", lw=1 )png = thispy.replace(".py", ".png")打印(写作",png )plt.savefig( png )plt.show()

I have some data such as:

I would like to fit a differentiable monotone curve to it. I tried PchipInterpolator class but on a very similar graph it resulted in:

This is not monotone.

How can I fit a monotone curve to data like this?

Here is a sample set of y values for another similar graph:

[0.1109157119023644, 0.20187393816931934, 0.14466318670239758, 0.16535159414166822, 0.05452708697483864, 0.2153046237959556, 0.2200300476272603, 0.21012762463269324, 0.15947100322395022, 0.2819691842129948, 0.15567770052985092, 0.24850595803020692, 0.1329341593280457, 0.15595107081606913, 0.3232021121832229, 0.23707961921686588, 0.2415887076540357, 0.32363506549779797, 0.3584089204036798, 0.29232772580068433, 0.22145994836140775, 0.22797587985241133, 0.2717787840603025, 0.3245255944762287, 0.29301098282789195, 0.32417076823344143, 0.3450906550996232, 0.34272097408024904, 0.3868714875012437, 0.41876692320045755, 0.3544198724867363, 0.33073960954801895, 0.3921033666371904, 0.33349050060172974, 0.3608862044547096, 0.37375822841635425, 0.5396399750708429, 0.4209201143798284, 0.42004773793166883, 0.5217725632679073, 0.5911731474218788, 0.43389609315065386, 0.4287288396176006, 0.43007525393257007, 0.5687062142675405, 0.6030811498722173, 0.5292225577714743, 0.47710974351051355, 0.6182720730381119, 0.6241033581931327, 0.6236788197617511, 0.6643161356364049, 0.5577616524049582, 0.6888440258481371, 0.6867893120660341, 0.6685257606057502, 0.599481675493677, 0.7309075091448749, 0.7644365338580481, 0.6176797601816733, 0.6751467827192018, 0.6452178017908761, 0.6684778262246701, 0.7003380077556168, 0.667035916425416, 0.8434451759113093, 0.8419343615815968, 0.8657695361433773, 0.7392487161484605, 0.8773282098364621, 0.8265679895117846, 0.7246599961191632, 0.7251899061730714, 0.9271640780410231, 0.9180581424305536, 0.8099033021701689, 0.8268585329594615, 0.8519967080830176, 0.8711231413093845, 0.8689802343798663, 0.8299523829217353, 1.0057741699770046, 0.8538130788729608, 0.9662784297225102, 1.023419780920539, 0.913146849759822, 0.9900885996579213, 0.8740638988529978, 0.8900285618419457, 0.9065474574434158, 1.0749522597307315, 1.0319120938258166, 1.0051369663172995, 0.9893558841613622, 1.051384986916457, 1.0327996870915341, 1.0945543972861898, 0.9716604944496021, 1.1490370559566179, 1.1379231481207432, 1.6836433783615088, 1.8162068766097395, 2.072155286917785, 2.0395966998366, 2.191064589600466, 2.1581974932543617, 2.163403843819597, 2.133441151300847, 2.1726053994136922, 2.1157865673629526, 2.2249636455682866, 2.2313062166802147, 2.1731708496472764, 2.315203950110816, 2.1601242661726827, 2.174940281421225, 2.2653635413275945, 2.337227057574145, 2.3645767548381618, 2.3084919291392527, 2.314014515926446, 2.25166717296155, 2.2621157708115778, 2.2644578546265586, 2.313504860292943, 2.398969190357051, 2.309443951779675, 2.278946047410807, 2.4080802287121146, 2.353652872018618, 2.35527529074088, 2.4233001060410784, 2.428767198055608, 2.35677123091093, 2.497135132404064, 2.3978099128437282, 2.3970802609341972, 2.4967434818740024, 2.511209192435555, 2.541001050440798, 2.5760248002036525, 2.5960512284192245, 2.4778408861721037, 2.5757724103530046, 2.631148267999664, 2.538327346218921, 2.4878734713248507, 2.6133797275761066, 2.6282561527857395, 2.6150327104952447, 3.102757164382848, 3.3318503012160905, 3.3907776288198193, 3.6065313558941936, 3.601180295875859, 3.560491539319038, 3.650095006265445, 3.574812155815713, 3.686227315374108, 3.6338261415040867, 3.5661194785086288, 3.5747332336054645, 3.560674343726918, 3.5678550481603635, 3.5342848534390967, 3.4929538312485913, 3.564544653619436, 3.6861775399566126, 3.6390300636595216, 3.6656336332413666, 3.5731185631923945, 3.5965520044069854, 3.537434489989021, 3.5590937423870144, 3.5331656424410083, 3.640652819618705, 3.5971240740252126, 3.641793843012055, 3.6064014089254295, 3.530378938786505, 3.613631139461306, 3.519542268056021, 3.5416251524576, 3.524789618934195, 3.5519951806099512, 3.6435695455293975, 3.6825670484650863, 3.5993379768209217, 3.628367553897596, 3.633290480934276, 3.5772841681579535, 3.602326323397947, 3.518180278272883, 3.531054006706696, 3.5566645495066167, 3.5410992153240985, 3.630762839301216, 3.5924649123201053, 3.646230633817883, 3.568290612034935, 3.638356129262967, 3.566083243271712, 3.6064978645771797, 3.4942864293427633, 3.595438454812999, 3.681726879126678, 3.6501308156903463, 3.5490717955938593, 3.598535359345363, 3.6328331698421654, 3.595159538698094, 3.556715819008055, 3.6292942886764554, 3.6362895697392856, 3.5965220100874093, 3.6103542985016266, 3.5715010140382493, 3.658769915445062, 3.5939686395400416, 3.4974461928859917, 3.5232691556732267, 3.6145687814416614, 3.5682054018341005, 3.648937250575395, 3.4912089018613384, 3.522426560340423, 3.6757968409374637, 3.651348691084845, 3.5395070091675973, 3.5306275536360383, 3.6153498246329883, 3.599762785949876, 3.5351931286962333, 3.6488316987683054, 3.5198301490992963, 3.5696570079786687, 3.561553836008927, 3.5659475947331423, 3.553147100256108, 3.5475591872743664, 3.6097226797553317, 3.6849600324757934, 3.5264731043844413, 3.506658609738451, 3.5535775980874114, 3.5487291053913554, 3.570651383823912, 3.552993371839188, 3.5054297764661846, 3.5723024888238792]
解决方案

Here's a monotone curve fitter in essentially 5 lines of python, with numpy and a lowpass filter from scipy.signal:

#!/usr/bin/env python
"""https://stackoverflow.com/questions/56551114/fully-monotone-interpolation-in-python """
# see also
# https://en.wikipedia.org/wiki/Monotone-spline aka I-spline
# https://scikit-learn.org/stable/modules/isotonic.html
# denis 2 March 2020

from __future__ import division, print_function
import numpy as np
from scipy import signal as sig

from matplotlib import pyplot as plt
import seaborn as sns


def butter_filtfilt( x, Wn=0.5, axis=0 ):
    """ butter( 2, Wn ), filtfilt
        axis 0 each col, -1 each row
    """
    b, a = sig.butter( N=2, Wn=Wn )
    return sig.filtfilt( b, a, x, axis=axis, method="gust" )  # twice, forward backward

def ints( x ):
    return x.round().astype(int)

def minavmax( x ):
    return "min av max %.3g %.3g %.3g" % (
            x.min(), x.mean(), x.max() )

def pvec( x ):
    n = len(x) // 25 * 25
    return "%s \n%s \n" % (
        minavmax( x ),
        ints( x[ - n : ]) .reshape( -1, 25 ))

#...............................................................................
def monofit( y, Wn=0.1 ):
    """ monotone-increasing curve fit """
    y = np.asarray(y).squeeze()
    print( "\n{ monofit: y %d %s  Wn %.3g " % (
        len(y), minavmax( y ), Wn ))
    ygrad = np.gradient( y )
    print( "grad y:", pvec( ygrad ))

        # lowpass filter --
    gradsmooth = butter_filtfilt( ygrad, Wn=Wn )
    print( "gradsmooth:", pvec( gradsmooth ))

    ge0 = np.fmax( gradsmooth, 0 )

    ymono = np.cumsum( ge0 )  # integrate, sensitive to first few
    ymono += (y - ymono).mean()

    err = y - ymono
    print( "y - ymono:", pvec( err ))
    errstr = "average |y - monofit|: %.2g" % np.abs( err ).mean()
    print( errstr )
    print( "} \n" )

    return ymono, err, errstr

#...............................................................................
if __name__ == "__main__":
    import sys

    np.set_printoptions( threshold=20, edgeitems=15, linewidth=120,
            formatter = dict( float = lambda x: "%.2g" % x ))  # float arrays %.2g
    print( 80 * "=" )

    thispy = sys.argv[0]
    infile = sys.argv[1] if len(sys.argv) > 1 \
        else "so-mono.txt"
    Wn = 0.1  # ?
    params = "%s %s  Wn %g " % (thispy, infile, Wn)
    print( params )

    y = np.loadtxt( infile ) * 100
    print( "y:", y )

    ymono, err, errstr = monofit( y, Wn=Wn )

    if 1:
        sns.set_style("whitegrid")
        fig, ax = plt.subplots( figsize=[10, 5] )
        plt.subplots_adjust( left=.05, right=.99, bottom=.05, top=.90 )
        fig.suptitle(
"Easy monotone curve fit: np.gradient | lowpass filter | clip < 0 | integrate \n"
            + errstr, multialignment="left" )

        ax.plot( ymono, color="orangered" )
        j = np.where( ymono < y )[0]
        xax = np.arange( len(y) )
        plt.vlines( xax[j], ymono[j], y[j], color="blue", lw=1 )
        j = np.where( ymono > y )[0]
        plt.vlines( xax[j], y[j], ymono[j], color="blue", lw=1 )

        png = thispy.replace( ".py", ".png" )
        print( "writing", png )
        plt.savefig( png )
        plt.show()

这篇关于python中的完全单调插值的文章就介绍到这了,希望我们推荐的答案对大家有所帮助,也希望大家多多支持!

09-23 01:23