zoukankan      html  css  js  c++  java
  • 一元二次曲线拟合Excel用

    import numpy
    import xlrd
    import matplotlib.pyplot as plt
    
    from scipy.optimize import leastsq
    
    X = numpy.array(range(0, 19))#2001,2019 
    # X = numpy.array(range(0, 30)) #x里面保存是0━29之间的所有整数
    
    # 读取excel数据
    def input_data(wbname):
        wb = xlrd.open_workbook(wbname)
        sh = wb.sheet_by_index(0)
                    
        table = [] # 存放每年被引频次
        sumy = [] # 总被引频次
        sum5 = [] # 发表最初五年被引频次
        #sum6_9 = [] # 随后四年???
        for row in range(1, sh.nrows):#逐行读取
            Y = sh.row_values(row)[32: 32 + 19]
            # Y = sh.row_values(row)[21: 21 + 30] #y里面保存的是这30年每年文献的被引次数。
            table.append(Y)
            sumy.append(sh.row_values(row)[19]) # 合计引用次数
            sum5.append(sh.row_values(row)[32: 32 + 5]) # 发表最初五年被引频次
            #sum6_9.append(sh.row_values(row)[32 + 5: 32 + 5 + 4]) # 随后四年
        return table, sumy, sum5
        #, sum6_9
    
    # 最小二乘法误差
    def RSquare(X, Y, f):
        yba = 0.0
        for y in Y:
            yba += y
        yba /= len(Y)
    
        SStot = 0.0
        for y in Y:
            SStot += (y - yba) ** 2
    
        SSres = 0.0
        for x, y in zip(X, Y):
            fx = f(x)
            SSres += (y - fx) ** 2
    
        return 1 - SSres / SStot
    
    if __name__ == '__main__':
        # table, sumlist = input_data("~/Desktop/毕业论文/引文/2001.xlsx")
        table, sumlist, sum5list = input_data("2001.xlsx")
    
        def residuals(p, x, y): # 定义计算误差函数
            polyfun = numpy.poly1d(p) # numpy.poly1d构造一个多项式
            return y - polyfun(x)
        
        def test_func(x, p1):
            a, b, c = p1
            return a*x**2+b*x+c
        
        #labels = []
        
        i = 2
        #j = 0
        for Y,sumy,sum5 in zip(table, sumlist,sum5list):
            #r = scipy.optimize.leastsq(residuals, [0, 0, 0], args=(X, numpy.array(Y))) # 最小二乘法
            r = leastsq(residuals, [0, 0, 0], args=(X, numpy.array(Y)))# 计算误差函数,拟合参数初始值,需要拟合的实验数据
            
            if r[1] not in [1, 2, 3, 4]:
                print("Error!")
                print(r)
                continue
            # print(r)
            p = r[0] # leastsq输出的第一个参数:拟合解组成的数组
            rsquare = RSquare(X, Y, numpy.poly1d(p))
    
            if rsquare >= 0.75 and sumy >= 20 and np.mean(sum5) <= 2: 
            # r方大于0.75,总被引频次大于20,发表最初五年年均被引频次不超过两次,随后四年至少为20次
                a, b, c = p
    
                t = -b / (2 * a)# 对称轴
                if 2.5 <= t <= 5.5: # 对称轴限定2.5 5.5
                    if a > 0:# 开口向上
                        print("{}:   {}   {}".format(i, rsquare, t)) # 第几行,r方,对称轴
                        print("{}:   {}   ".format(i, p)) # 第几行,拟合的参数
                        #print(numpy.poly1d(p)) # numpy.poly1d(p)是一个函数
                                          
                        colorlist=plt.cm.cool(np.linspace(0,1,7))
                        markers=["D","x","h",".","^",">","v"]
                        #labels.append(r'$y =%ix**2 + %ix + %i$'% (p[0],p[1],p[2]))
                        
                        for n in range(0,4):# 几张图
                            plt.plot(X, numpy.array(Y),color=colorlist[2],marker=markers[0],linestyle="-") # 被引频次图,蓝色是真实的,橙色的是拟合出来的
                            plt.plot(X, test_func(X, p),color=colorlist[2],marker=markers[2],linestyle="--") # 一元二次函数
                            #plt.legend(labels)
                        plt.show()
                        #j+=1 # 改变形状颜色
                        
            i += 1
    

    代码部分参考学妹毕业论文。

  • 相关阅读:
    MySQL 5.6.9 RC 发布
    红薯 Java 8 的日期时间新用法
    Couchbase Server 2.0 发布,NoSQL 数据库
    Firefox OS 模拟器 1.0 发布
    Calculate Linux 13 Beta 1 发布
    敏捷测试的团队构成
    Node.js 0.8.16 发布(稳定版)
    JASocket 1.1.0 发布
    Samba 4.0 正式版发布,支持活动目录
    Seafile 1.3 发布,文件同步和协作平台
  • 原文地址:https://www.cnblogs.com/Cookie-Jing/p/15427517.html
Copyright © 2011-2022 走看看