python - 评估 numpy 多项式返回错误的值

Question

我有一个系数为 c=2.44717404e-03 和 m=1.88697661e+01 的线性 numpy 多项式，但是当我尝试对其求值时，它给出了错误的输出。

f=np.polynomial.polynomial.Polynomial([2.44717404e-03, 1.88697661e+01])

例如

f(0.015)

返回

3.9646796502044523

而不是正确的结果

0.28549366554

# -*- coding: utf-8 -*-
"""
Created on Fri Mar 13 17:15:51 2020

@author: tomhe
"""
#Modules
if True:
    import numpy as np
    import matplotlib.pyplot as plt
    import pandas as pd
    from scipy import stats
    import scipy.optimize as scpo

#Functions
def M(params, *args):
    c=params[0]
    m=params[1]

    x=args[0]
    y=args[1]
    ex=args[2]
    ey=args[3]

    return(sum((y-(m*x+c))**2/(ey**2+m**2*ex**2)))

#Data
if True:
    df=pd.read_excel("Biot Savart Experiment Data 2.xlsx")

    B=np.array(df['B'].values)
    x=np.arange(0,18,0.5)
    a=3.95
    x1star=6.35
    x2star=10.3
    I=20.02

    eB=np.array([0.001]*len(B))
    ex=0.05
    ea=0.1
    ex1star=0.05
    ex2star=0.05
    eI=0.01

#Linearisation
if True:
    z=((x-x1star)**2+a**2)**(-3/2)+((x-x2star)**2+a**2)**(-3/2)
    ez=3*np.sqrt(((x-x1star)**2+a**2)**(-5)*((x-x1star)**2*(ex**2+ex1star**2)+(a*ea)**2)+((x-x2star)**2+a**2)**(-5)*((x-x2star)**2*(ex**2+ex2star**2)+(a*ea)**2))

#Linear regression
if True:
    optrslt=scpo.minimize(M,[0,20],(z,B,ez,eB))
    optparams=optrslt.x
    f=np.polynomial.polynomial.Polynomial(optparams,domain=(min(z),max(z)))
    linfitcoords=f.linspace(1000)

#Resiudals
if True:
    residuals=f(z)-B

#Residual level line
if True:
    rf=np.polynomial.polynomial.Polynomial([0],domain=(min(z),max(z)))
    rfcoords=rf.linspace(1000)

#Plotting
if True:
    fig=plt.figure(1)
    frame=plt.gca()
    axis1=fig.add_axes((0,0,1,1))
    axis1.errorbar(z,B,xerr=ez,fmt='o', ms=1, lw=0.5, color='red', capsize=2, capthick=0.5, ecolor='black')
    #axis1.plot(z,line(z,optparams[1],optparams[0]),lw=0.5,color='blue')
    #axis1.plot(x,y,lw=0.5,color='blue')
    axis1.plot(linfitcoords[0],linfitcoords[1],lw=0.5,color='blue')
    plt.ylabel('Magnetic Field Strength ($mT$)')
    axis1.xaxis.set_visible(False)
    axis1.set_xlim(left=0,right=0.025)
    axis1.set_ylim(bottom=0,top=0.5)
    axis2=fig.add_axes((0,-0.4,1,0.4))
    axis2.set_xlim(left=0,right=0.025)
    axis2.plot(rfcoords[0],rfcoords[1],lw=0.5, color='blue')
    axis2.scatter(z,residuals, s=1, color='red')
    plt.ylabel("Residuals ($mT)$")
    axis1.yaxis.set_label_coords(-0.1,0.5)
    axis2.yaxis.set_label_coords(-0.1,0.5)
    plt.xlabel('z $(cm^{-3})$')
    plt.xlim([0,0.025])
    axis3=fig.add_axes((1,-0.4,0.2,0.4))
    axis3.hist(residuals,orientation='horizontal',color='gray')
    axis3.yaxis.set_visible(False)
    axis3.xaxis.tick_top()
    axis3.xaxis.set_label_position('top') 
    plt.xlabel('Frequency')
    plt.xticks([1,3,5,7])
    #Note:restarting kernel (close window) fixes x axis label problem where you accidentally set as variable and then it breaks

#Save Figure
if True:
    plt.savefig("Linearised_Biot_Savart_Graph.png",dpi=1000,bbox_inches='tight')

#Gradient and mu_0 analysis
if True:
    m=optparams[1]
    mu0=2*m/(I*a**2)
    em=1#need to actually figure this out
    emu0=mu0*np.sqrt((em/m)**2+(eI/I)**2+4*(ea/a)**2)
   # print(mu0)
    #print(emu0)

#Normality of residuals
if True:
    print(stats.shapiro(residuals))
    print(stats.normaltest(residuals))
    #print(stats.kstest(residuals,stats.norm(np.mean(residuals),np.var(residuals)).cdf))

注意:所有的 if 语句只是为了方便 Spyder 中的代码折叠。

就上下文而言，我正在对磁场的一些物理数据进行数据分析和可视化。我的主要问题是我试图绘制一条最适合数据的线，即使我得到的参数对于最适合的线是正确的，当我尝试使用 numpy 多项式生成坐标时绘制最佳拟合线，它会生成错误的坐标。

Answer 1

看起来问题是传递给 f 的 domain 参数

来自docs您正在传递域的值但没有传递窗口的值，因此 f 根据 map [min(z), max(z)] -> [-1, 1] 映射输入]

这就是输出不同的原因。

要清楚地看到这一点，请看下面的例子

import numpy as np

# x**2 with no domain shift
f = np.polynomial.polynomial.Polynomial([0,0,1])
print(f)
>> poly([0. 0. 1.])
print(f(1))
>> 1.0

# x**2 with domain shift [-0.1, 0.1] -> [-1, 1], so an input of 1 maps to 10
g = np.polynomial.polynomial.Polynomial([0,0,1], domain=[-0.1, 0.1])
print(g)
>> poly([0. 0. 1.])
print(g(1))
>> 100.0

因此任何不是 [-1,1] 的 domain 值都会给出一个多项式，由于输入的缩放而返回“错误”的输出