一元二次曲線擬合Excel用

對excel進行操作，有別於csv，具體的表現在excel是一個工作簿，csv只有一個表。詳細的需要進一步整理。

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

代碼部分參考學妹畢業論文。