穿

穿山aA

V1

2023/04/10阅读:35主题:科技蓝

多项式拟合和龚珀兹曲线

多项式预测和龚珀兹曲线

import numpy as np
from scipy.stats import norm
from scipy.stats import t
import matplotlib.pyplot as plt
#!/usr/bin/python
# -*- encoding: utf-8 -*-
'''
@File    :   多项式模型的建立.ipynb
@Time    :   2023/04/10 13:59:05
@Author  :   chuanshana
@Contact :   算了别联系
'''

import matplotlib.pyplot as plt 
font = {'family' : 'KaiTi','weight' : 'bold','size'   : '16'}
plt.rc('font', **font)               # 步骤一(设置字体的更多属性)
plt.rc('axes', unicode_minus=False)  # 步骤二(解决坐标轴负数的负号显示问题)
# 设置图形显示风格
plt.style.use('ggplot')

先定义一个求差分的function

def get_difference(x,a):
    b=1
    if  a==1:
        difference=[x[i]-x[i-1]for i in range(1,len(x))]
        return  ["----"]+difference
    df=[x[i]-x[i-1]for i in range(1,len(x))]
    while b<a:
        answer=[df[i]-df[i-1]for i in range(1,len(df))]
        df=answer
        b+=1
    return ["----"]*(len(x)-len(df))+df
y=[2790,2950,3140,3350,3588,3862,4168]
get_difference(y,1)   

​ ['----', 160, 190, 210, 238, 274, 306]

第一题多项式预测

画个图看看是什么曲线,在求差分看看

import pandas as pd
year=[i for i in range(2008,2015)]
y_t=[2790,2950,3140,3350,3588,3862,4168]
#首先画个图看看啥情况
plt.plot(year,y)#这个图看着像一次的,进行差分看看
差分计算表=pd.DataFrame({
    "y_t":y_t,
    "一阶差分":get_difference(y,1),
    "二阶差分":get_difference(y,2),
    "三阶差分":get_difference(y,3),
})
差分计算表
一阶差分 二阶差分 三阶差分
0 2790 ---- ---- ----
1 2950 160 ---- ----
2 3140 190 30 ----
3 3350 210 20 -10
4 3588 238 28 8
5 3862 274 36 8
6 4168 306 32 -4
  • 从图片看像是线性,但是眼睛看的不一定真实,进行差分运算,判断是二次

计算二次多项式曲线模型参数

timeseries=[i for i in range(-3,4)]
t_square=[i**2 for i in timeseries]
t_4=[i**4 for i in timeseries];t_4
ty=[i*j for i,j in zip(timeseries,y_t)]
t_square_y=[i**2*j for i,j in zip(timeseries,y_t)]
dataframe=pd.DataFrame({"年份":year,
                        "timeseries":timeseries,
                       "y_t":y_t,
                       "t_square":t_square,
                       "t_4":t_4,
                       "ty":ty,
                       "t_square_y":t_square_y}
                       )
print("下面是二次多项式曲线模型参数计算表")
dataframe=dataframe.set_index("年份")
dataframe.loc['Row_sum'] = dataframe.apply(lambda x: x.sum())
dataframe

​ [81, 16, 1, 0, 1, 16, 81]

​ 下面是二次多项式曲线模型参数计算表

时序
年份
2008 -3 2790 9 81 -8370 25110
2009 -2 2950 4 16 -5900 11800
2010 -1 3140 1 1 -3140 3140
2011 0 3350 0 0 0 0
2012 1 3588 1 1 3588 3588
2013 2 3862 4 16 7724 15448
2014 3 4168 9 81 12504 37512
Row_sum 0 23848 28 196 6406 96598
  • 带入三元一次方程组 很容易解得:[{b0: 23446/7, b2: 201/14, b1: 3203/14}] 这sb东西把b1放后面干什么,抄错了沃日

#解方程
from sympy import *
b0,b1,b2=symbols("b0,b1,b2")
eq1=7*b0+28*b2-23848
eq2=28*b1-6406
eq3=28*b0+196*b2-96598
solve([eq1,eq2,eq3],[b0,b1,b2],dict=True)
b0=23446/7
b1=3203/14
b2=201/14

​ [{b0: 23446/7, b2: 201/14, b1: 3203/14}]

  • 若要预测2105年的捕捞量,则 时:

eee,跟标准答案差了140多,不过由误差正常

估计误差计算表

y_hat=[round(23446/7+3203/14*i+201/14*i**2,3)for i in timeseries]
e_t=[i-j for i,j in zip(y,y_hat)]
e_t_square=[(i-j)**2 for i,j in zip(y,y_hat)]
估计误差计算表=pd.DataFrame({
    "y":y,
    "年份":year,
    "y_hat":y_hat,
    "(y-y_hat)":e_t,
    "(y-y_hat)^2":e_t_square
}).set_index("年份")
估计误差计算表
y
年份
2008 2790 2792.286 -2.286 5.225796
2009 2950 2949.286 0.714 0.509796
2010 3140 3135.000 5.000 25.000000
2011 3350 3349.429 0.571 0.326041
2012 3588 3592.571 -4.571 20.894041
2013 3862 3864.429 -2.429 5.900041
2014 4168 4165.000 3.000 9.000000
SE=sqrt(估计误差计算表["(y-y_hat)^2"].sum()/(7-3));SE
y_2015=b0+b1*4+b2*16;y_2015
置信区间=[y_2015-t.ppf(0.975,6)*SE,y_2015+t.ppf(0.975,6)*SE]
置信区间

原来是这个傻逼答案的问题,我写的没问题

第二题龚珀(po第四声)兹曲线

  • 这狗几把不是取以e为底的对数,而是以10为底,即 我真服了

  • 第一步还是画个图,看起来长得真像线性我服了

year=[i for i in range(2009,2015)]
销售量=[8.7,10.6,13.3,16.5,20.6,26]
plt.plot(year,销售量)
#定义一个函数,往pandas中任意一行添加数据
def insert(df, i, df_add):
    # 指定第i行插入一行数据
    df1 = df.iloc[:i, :]
    df2 = df.iloc[i:, :]
    df_new = pd.concat([df1, df_add, df2], ignore_index=True)
    return df_new
lgy=np.log10(销售量)
销售量情况=pd.DataFrame({
    "年份":year,
    "时序(t)":[i for i in range(0,6)],
     "销售量(y)":销售量,
    "lgy":lgy
}).set_index("年份")

2个lgy求和=sum(lgy[0:2])
中间2个lgy求和=sum(lgy[2:4])
2个lgy求和=sum(lgy[4:6])

sumIlgy=["--","--",1]
df_add1 = pd.DataFrame({'lgy':[前2个lgy求和]})
df_add2 = pd.DataFrame({'lgy':[中间2个lgy求和]})
#在第2行插入一条新的数据
销售量情况= insert(销售量情况, 2, df_add1)
销售量情况= insert(销售量情况, 5, df_add2)
#最后一行添加一个最后两个求和
销售量情况.loc["disange"]={
    "时序(t)":----,
     "销售量(y)":----,
    "lgy":后2个lgy求和}
index=[2009,2010,"sumIlgy",2011,2012,"sumIIlgy",2013,2014,"sumIIIlgy"]
销售量情况["年份"]=index
销售量情况.set_index("年份")
时序(t) 销售量(y)
年份
2009 0.0 8.7 0.939519
2010 1.0 10.6 1.025306
sumIlgy NaN NaN 1.964825
2011 2.0 13.3 1.123852
2012 3.0 16.5 1.217484
sumIIlgy NaN NaN 2.341336
2013 4.0 20.6 1.313867
2014 5.0 26.0 1.414973
sumIIIlgy NaN NaN 2.728841

分类:

后端

标签:

后端

作者介绍

穿
穿山aA
V1