4

我的数据点很少,我想创建一条线以在以符号学比例绘制时最适合数据点。我已经尝试过 scipy 的曲线拟合和三次插值,但与数据趋势相比,它们对我来说似乎都不是很合理。

我恳请您检查是否有更有效的方法来创建适合数据的直线。可能外推可以做,但我没有找到关于python外推的好文档。

非常感谢您的帮助

import sys
import os
import numpy
import matplotlib.pyplot as plt
from pylab import *
from scipy.optimize import curve_fit
import scipy.optimize as optimization
from scipy.interpolate import interp1d
from scipy import interpolate

在此处输入图像描述

Mass500 = numpy.array([ 13.938 , 13.816,  13.661,  13.683,  13.621,  13.547,  13.477,  13.492,  13.237,
  13.232,  13.07,   13.048,  12.945,  12.861,  12.827,  12.577,  12.518])

y500 = numpy.array([  7.65103978e-06,   4.79865790e-06,   2.08218909e-05,   4.98385924e-06,
   5.63462673e-06,   2.90785458e-06,   2.21166794e-05,   1.34501705e-06,
   6.26021870e-07,   6.62368879e-07,   6.46735547e-07,   3.68589447e-07,
   3.86209019e-07,   5.61293275e-07,   2.41428755e-07,   9.62491134e-08,
   2.36892162e-07])



plt.semilogy(Mass500, y500, 'o')


# interpolation
f2 = interp1d(Mass500, y500, kind='cubic')
plt.semilogy(Mass500, f2(Mass500), '--')


# curve-fit
def line(x, a, b):
    return 10**(a*x+b)

#Initial guess.
x0     = numpy.array([1.e-6, 1.e-6])

print optimization.curve_fit(line, Mass500, y500, x0)
popt, pcov = curve_fit(line, Mass500, y500)
print popt
plt.semilogy(Mass500, line(Mass500, popt[0], popt[1]), 'r-')



plt.legend(['data', 'cubic', 'curve-fit'], loc='best')

show()
4

2 回答 2

5

numpy 和 scipy 中有许多回归函数可用。 scipy.stats.lingress是更简单的函数之一,它返回常见的线性回归参数。

以下是拟合半对数数据的两个选项:

  1. 绘制转换后的数据
  2. 重新缩放轴并转换输入/输出函数值

给定

import numpy as np
import scipy as sp
import matplotlib.pyplot as plt

%matplotlib inline


# Data
mass500 = np.array([ 
    13.938 , 13.816,  13.661,  13.683,
    13.621,  13.547,  13.477,  13.492,
    13.237,  13.232,  13.07,   13.048,  
    12.945,  12.861,  12.827,  12.577,  
    12.518
])

y500 = np.array([  
    7.65103978e-06,   4.79865790e-06,   2.08218909e-05,   4.98385924e-06,
    5.63462673e-06,   2.90785458e-06,   2.21166794e-05,   1.34501705e-06,
    6.26021870e-07,   6.62368879e-07,   6.46735547e-07,   3.68589447e-07,
    3.86209019e-07,   5.61293275e-07,   2.41428755e-07,   9.62491134e-08,
    2.36892162e-07
])

代码

选项 1:绘制转换后的数据

# Regression Function
def regress(x, y):
    """Return a tuple of predicted y values and parameters for linear regression."""
    p = sp.stats.linregress(x, y)
    b1, b0, r, p_val, stderr = p
    y_pred = sp.polyval([b1, b0], x)
    return y_pred, p

# Plotting
x, y = mass500, np.log(y500)                      # transformed data
y_pred, _ = regress(x, y)

plt.plot(x, y, "mo", label="Data")
plt.plot(x, y_pred, "k--", label="Pred.")
plt.xlabel("Mass500")
plt.ylabel("log y500")                            # label axis
plt.legend()

输出

在此处输入图像描述

一种简单的方法是绘制转换后的数据并标记适当的对数轴。


选项 2:重新缩放轴并转换输入/输出函数值

代码

x, y = mass500, y500                              # data, non-transformed
y_pred, _ = regress(x, np.log(y))                 # transformed input             

plt.plot(x, y, "o", label="Data")
plt.plot(x, np.exp(y_pred), "k--", label="Pred.") # transformed output
plt.xlabel("Mass500")
plt.ylabel("y500")
plt.semilogy()
plt.legend()

输出

在此处输入图像描述

第二种选择是将轴更改为半对数刻度(通过plt.semilogy())。在这里,未转换的数据自然呈现线性。另请注意,标签按原样表示数据。

为了进行准确的回归,剩下的就是将传递给回归函数的数据(通过np.log(x)np.log10(x))转换为返回正确的回归参数。当使用互补运算(即np.exp(x)或)绘制预测值时,这种转换会立即反转10**x

于 2018-01-28T03:36:51.690 回答
3

如果您想要一条在 log-y 尺度上看起来不错的线,则将一条线拟合到 y 值的对数。

def line(x, a, b):
    return a*x+b
popt, pcov = curve_fit(line, Mass500, np.log10(y500))
plt.semilogy(Mass500, 10**line(Mass500, popt[0], popt[1]), 'r-')

就是这个; 我只省略了看起来不相关的三次插值部分。

合身

于 2018-01-28T03:43:22.673 回答