沉降物化實驗的python數據分析與可視化

2021-02-21 波本的咖啡屋
1. 引言

上期做無機實驗的數據分析可視化的時候預告過這一期，但前段時間肝各路論文，也沒抽出太多時間來寫，今天補上。

沉降物化實驗是通過測定一個多級分散體系沉降曲線，得到分散體系粒度分布的實驗，其中會要求使用Origin進行數據處理，因為這個實驗涉及到幾個比較麻煩的地方

根據實驗所給出的沉降曲線數據(m-t圖)，繪製沉降曲線並根據已知超越方程對曲線進行擬合通過擬合所得的方程，求導函數，求取給定粒徑的粒子沉降完全所對應的時間點並通過做切線的方法得到此時給定粒徑區間的粒子的沉降量

用python進行曲線擬合是我此前一直沒解決的問題，此前我比較喜歡打開matlab的curvefit工具箱進行擬合，不得不說嘗試性的曲線擬合確實是用那個方法比較好，但通過在CSDN學習，我發現scipy.optimize模塊下有個curve_fit對象，可以用來進行曲線擬合，配合spyder可以達到與matlab的curvefit工具包無異的效果

2. 具體數據處理方法2.1 實驗數據導入

首先導入實驗數據，目前還沒有做成項目的樣子，所以導入方法很直接。csv文件是根據實驗導出的數據文件經過一次簡單加工而成，是進行了零點較正的數據。

# 引入所需庫
import numpy as np
import matplotlib.pyplot as plt
import pandas as pd
import statsmodels.api as sm
from statistic_model import linear  # 引入線性擬合
from scipy.optimize import curve_fit

# 數據導入
data = pd.read_csv('10183791.csv')
y = np.array(data['m'])
x = np.array(data['t'])
# 實驗參數
eta = 0.0009111  # 黏度,Pa·S
rho = 2700 # 滑石粉密度,kg/m3
rho_0 = 997.2995 # 介質密度,kg/m3
h = 0.13 # 沉降高度,m
g = 9.81
至於為什麼要引入線性擬合模塊？還請客官往下看
2.2 沉降曲線的擬合2.2.1 沉降曲線的表達形式
根據實驗講義，沉降曲線的擬合方程為：
此處一個大問題是
取最後八組數據，求1/t。將m對1/t線性擬合，得到一條擬合直線，這條擬合直線的R^2應當是0.95以上的線性擬合外推到1/t=0，對應的m即為曲線上界

# 線性擬合求mc
x_c = x[-8:] # t
y_c = y[-8:] # m
x_cc = 1000/x_c
P_linear, R2 = linear(x_cc, y_c)
xx_c = np.linspace(np.min(x_cc), np.max(x_cc), 50)
yy_c = np.polyval(P_linear, xx_c)

mc = P_linear[1] #最大沉降質量,g
print(f'mc:{mc}')

# 線性擬合圖作圖
# 作圖
plt.figure(1)
plt.rcParams['font.family'] = 'SimHei'
plt.rcParams['font.size'] = 16
# 調整matplotlib內部設置使之支持中文,這種方法改變全局設置
plt.plot(xx_c, yy_c, 'r-.', label='擬合曲線', linewidth=4)
plt.plot(x_cc, y_c, 'bo', label='實驗數據點', markersize=10)
plt.xlabel(u'$1000/t (s^{-1})$', fontsize=24)
plt.ylabel(u'沉降重量 $m/g$', fontsize=24)
plt.title('線性擬合外推法求取沉降總量曲線圖', fontsize=28)
plt.grid()
plt.legend()
# 擬合線描述圖例
txt_point = (np.mean(xx_c[xx_c >= np.mean(xx_c)]), np.mean(yy_c))
# 相對點法
msg = '擬合方程: y={:.4g}x+{:.4g}\n  R2 = {:.4f}'.format(
    P_linear[0], P_linear[1], R2)
plt.annotate(msg, xy=txt_point, xytext=txt_point)
即可以通過線性擬合的方法求得
由於之前的Arial Unicode MS在書寫"沉"字的時候看起來很醜，因此將字體換成黑體SimHei
線性擬合曲線圖
接下來是重頭戲，基於已知曲線進行曲線擬合
# 曲線擬合函數定義
def func(t,a,b,c):
    global mc
    return mc * (1-np.exp(-a * t**(b + c*np.log(t))))

# curve_fit環節
popt,pcov = curve_fit(func,x,y)
a,b,c = popt[0],popt[1],popt[2]
yy = func(x,a,b,c)
print(a,b,c)

# 作圖
plt.figure(2)
plt.rcParams['font.family'] = 'SimHei'
plt.rcParams['font.size'] = 16
# 調整matplotlib內部設置使之支持中文,這種方法改變全局設置
plt.plot(x, yy, 'r-', label='擬合曲線', linewidth=2.2)
plt.plot(x, y, 'bo', label='實驗數據點', markersize=5)
plt.xlabel('時間 $t/s$', fontsize=24)
plt.ylabel('沉降重量 $m/g$', fontsize=24)
plt.title('多級分散體系沉降曲線圖', fontsize=28)
# 擬合線描述圖例
txt_point = (np.mean(x), np.mean(yy))
# 相對點法
msg = '擬合方程:\n $m_t = mc * [1-exp(- (a) * t ^{(b + c \ln (t)) } ) ]$'
msg = msg.replace('mc', f'{mc:.4g}')
msg = msg.replace('a', f'{a:.4g}')
msg = msg.replace('b', f'{b:.4g}')
msg = msg.replace('c', f'{c:.4g}')
plt.annotate(msg, xy=txt_point, xytext=txt_point)
plt.grid()
plt.legend()
得到曲線擬合結果和沉降曲線圖：
沉降曲線圖2.2.2 基於擬合結果的進一步數據處理
在使用scipy.optimize.curve_fit的時候，將函數單獨使用def定義出來用，這樣看起來代碼多了，實際上是將函數定義和函數使用分離，算是一種模塊化的思想。接下來我們也打算用這種思路進行數據處理：
# 導函數
# 導函數形式是固定的 直接內置
def dfunc(t,a,b,c):
    global mc
    return a*mc*t**(b + c*np.log(t))*(
        (b + 2*c*np.log(t))/t)*np.exp(-a*t**(b + c*np.log(t)))

# 由半徑求取對應的完全沉降時間的函數
def func_tr(r, eta=eta, rho=rho, rho_0=rho_0, h=h, g=g):
    t = (9*eta*h)/(2*g*(rho-rho_0)*r**2)
    return t

# 利用切線斜率求截距
def func_St(t, mt, dmdt):
    S = mt - t*dmdt
    return S
# 數據處理求分布函數折線圖
r = np.arange(1,14) # um
dr = 1 # um
t_r = func_tr(r*1E-6)
m_tr = func(t_r,a,b,c)
dmdt_r = dfunc(t_r,a,b,c)
S_r = func_St(t_r, m_tr, dmdt_r)
dS_r = []
for i in range(1,len(S_r)):
    dS_r.append(S_r[i-dr]-S_r[i])
dS_r = np.array(dS_r)
f_Sr = dS_r / (mc * dr * 1E-6)
可以得到一系列的輸出數據，這些輸出數據可以在spyder裡面直接可視化，我在代碼末尾設計了一段輸出，用來輸出這些數據處理的所得結果。這裡的數據處理過程可以參照沉降物化的實驗講義
# 出數據
result_data = pd.DataFrame(
    [r,t_r,m_tr,dmdt_r,S_r,dS_r,f_Sr],
    index=['r','t','m_t','dm/dt','S','∆S','f'])
result_data.to_csv('result_cjwh.csv')
# 出圖
plt.show()

2.2.3 粒度分布圖輸出
# 作圖
plt.figure(3)
plt.rcParams['font.family'] = 'SimHei'
plt.rcParams['font.size'] = 16
# 調整matplotlib內部設置使之支持中文,這種方法改變全局設置
plt.bar(r[1:], f_Sr, width=1,edgecolor='red',
    color='white',linewidth=4)
plt.xlabel('半徑 $r / um$', fontsize=24)
plt.ylabel('分布函數 f', fontsize=24)
plt.title('多級沉降體系粒度分布圖', fontsize=28)
基本的柱形圖做法即可，關於plt.bar可以在CSDN上進一步了解
粒度分布圖3. 全部原始碼
# JamesBourbon in 20201116
# 沉降物化實驗的數據處理
# 曲線擬合scipy.optimize.curve_fit初應用

import numpy as np
import matplotlib.pyplot as plt
import pandas as pd
import statsmodels.api as sm
from statistic_model import linear  # 引入線性擬合
from scipy.optimize import curve_fit

# 數據導入
data = pd.read_csv('10183791.csv')
y = np.array(data['m'])
x = np.array(data['t'])
# 實驗參數
eta = 0.0009111  # 黏度,Pa·S
rho = 2700 # 滑石粉密度,kg/m3
rho_0 = 997.2995 # 介質密度,kg/m3
h = 0.13 # 沉降高度,m
g = 9.81

# 線性擬合求mc
x_c = x[-8:] # t
y_c = y[-8:] # m
x_cc = 1000/x_c
P_linear, R2 = linear(x_cc, y_c)
xx_c = np.linspace(np.min(x_cc), np.max(x_cc), 50)
yy_c = np.polyval(P_linear, xx_c)

mc = P_linear[1] #最大沉降質量,g
print(f'mc:{mc}')

# 線性擬合圖作圖
# 作圖
plt.figure(1)
plt.rcParams['font.family'] = 'SimHei'
plt.rcParams['font.size'] = 16
# 調整matplotlib內部設置使之支持中文,這種方法改變全局設置
plt.plot(xx_c, yy_c, 'r-.', label='擬合曲線', linewidth=4)
plt.plot(x_cc, y_c, 'bo', label='實驗數據點', markersize=10)
plt.xlabel(u'$1000/t (s^{-1})$', fontsize=24)
plt.ylabel(u'沉降重量 $m/g$', fontsize=24)
plt.title('線性擬合外推法求取沉降總量曲線圖', fontsize=28)
plt.grid()
plt.legend()
# 擬合線描述圖例
txt_point = (np.mean(xx_c[xx_c >= np.mean(xx_c)]), np.mean(yy_c))
# 相對點法
msg = '擬合方程: y={:.4g}x+{:.4g}\n  R2 = {:.4f}'.format(
    P_linear[0], P_linear[1], R2)
plt.annotate(msg, xy=txt_point, xytext=txt_point)

# 曲線擬合函數定義
def func(t,a,b,c):
    global mc
    return mc * (1-np.exp(-a * t**(b + c*np.log(t))))

# curve_fit環節
popt,pcov = curve_fit(func,x,y)
a,b,c = popt[0],popt[1],popt[2]
yy = func(x,a,b,c)
print(a,b,c)

# 作圖
plt.figure(2)
plt.rcParams['font.family'] = 'SimHei'
plt.rcParams['font.size'] = 16
# 調整matplotlib內部設置使之支持中文,這種方法改變全局設置
plt.plot(x, yy, 'r-', label='擬合曲線', linewidth=2.2)
plt.plot(x, y, 'bo', label='實驗數據點', markersize=5)
plt.xlabel('時間 $t/s$', fontsize=24)
plt.ylabel('沉降重量 $m/g$', fontsize=24)
plt.title('多級分散體系沉降曲線圖', fontsize=28)
# 擬合線描述圖例
txt_point = (np.mean(x), np.mean(yy))
# 相對點法
msg = '擬合方程:\n $m_t = mc * [1-exp(- (a) * t ^{(b + c \ln (t)) } ) ]$'
msg = msg.replace('mc', f'{mc:.4g}')
msg = msg.replace('a', f'{a:.4g}')
msg = msg.replace('b', f'{b:.4g}')
msg = msg.replace('c', f'{c:.4g}')
plt.annotate(msg, xy=txt_point, xytext=txt_point)
plt.grid()
plt.legend()

# 導函數
# 導函數形式是固定的 直接內置
def dfunc(t,a,b,c):
    global mc
    return a*mc*t**(b + c*np.log(t))*(
        (b + 2*c*np.log(t))/t)*np.exp(-a*t**(b + c*np.log(t)))

# 由半徑求取對應的完全沉降時間的函數
def func_tr(r, eta=eta, rho=rho, rho_0=rho_0, h=h, g=g):
    t = (9*eta*h)/(2*g*(rho-rho_0)*r**2)
    return t

# 利用切線斜率求截距
def func_St(t, mt, dmdt):
    S = mt - t*dmdt
    return S
# 數據處理求分布函數折線圖
r = np.arange(1,14) # um
dr = 1 # um
t_r = func_tr(r*1E-6)
m_tr = func(t_r,a,b,c)
dmdt_r = dfunc(t_r,a,b,c)
S_r = func_St(t_r, m_tr, dmdt_r)
dS_r = []
for i in range(1,len(S_r)):
    dS_r.append(S_r[i-dr]-S_r[i])
dS_r = np.array(dS_r)
f_Sr = dS_r / (mc * dr * 1E-6)

# 作圖
plt.figure(3)
plt.rcParams['font.family'] = 'SimHei'
plt.rcParams['font.size'] = 16
# 調整matplotlib內部設置使之支持中文,這種方法改變全局設置
plt.bar(r[1:], f_Sr, width=1,edgecolor='red',
    color='white',linewidth=4)
plt.xlabel('半徑 $r / um$', fontsize=24)
plt.ylabel('分布函數 f', fontsize=24)
plt.title('多級沉降體系粒度分布圖', fontsize=28)

# 出數據
result_data = pd.DataFrame(
    [r,t_r,m_tr,dmdt_r,S_r,dS_r,f_Sr],
    index=['r','t','m_t','dm/dt','S','∆S','f'])
result_data.to_csv('result_cjwh.csv')
# 出圖
plt.show()

4. 總結
這個實驗的實驗結果數據文件少則幾十行，多的可以上千行，但都可以囊括在這套代碼裡面，出圖的時候調一下樣品點和曲線的大小就行。這套代碼幾乎把實驗級別的數據分析可視化做到極致了，再往上就要大數據和資料庫人工智慧級別了。
在這邊講解的時候沒有涉及到具體的各個庫各個方法各各個函數的使用細節，只是說了需要解決的問題和解決問題的方法，如果大家有這方面的需求歡迎後臺留言，考慮一下2333
事情好多，昨天才忙完成思危名譽校長獎學金答辯和創新實驗大賽論文，今天又要寫實驗報告忙論文改代碼跑計算，總歸事情是做不完的，還在更公眾號真的是熱情了，B站和雲峰聯的事情都只能寒假再說的樣子😂
不過我會繼續努力的吧
希望真的能如自己在成獎答辯的時候說的一樣，這一年能做到
鯤鵬展翅，直上九天。
歡迎點讚收藏在看三連呀，多多指教~

相關焦點

  • 【Python教程】用Python進行數據可視化
    ，我們可以通過數據分析信息是可量化的、可測度的、有形式的，可被訪問、生成、存儲、分發、搜索、壓縮和複製。當數據發生變化時，這種演變過程隨之發生。通過離散算法　     數據的種類不同，包括表現數據, 實驗數據,基準數據2.2 可視化的思路和工具可視化的整個過程需要不用技能和專業領域的人。
  • python數據分析專題 (7):python數據分析模塊
    python是一門優秀的程式語言，而是python成為數據分析軟體的是因為python強大的擴展模塊。
  • 13分鐘,教你python可視化分析20W數據,找到妹子最愛的內衣
    我打開某寶首頁，頓時想對某些商品的評論進行一次爬取分析。一提到爬取，我們肯定先想到python，那可真的是利器，不過我提醒一句，不允許爬的千萬不要碰。python雖然爬取簡單，但是做可視化分析需要費不少力氣，雖然有 Matplotlib 和 Seaborn 兩個包就足夠了。
  • Python數據分析可視化教程學習
    來源：博學谷 作者：照照隨著數據分析與數據挖掘行業的日益火爆，Python數據分析可視化教程也成為了許多人提升職業競爭力的必修課。本文將詳細為大家介紹Python數據分析可視化教程的具體內容，全文會圍繞著為什麼學？學什麼？學了有什麼用？
  • 推薦: 一本「高顏值」的Python語言數據可視化圖書
    在先介紹數據可視化基礎理論後，本書系統性地介紹了幾乎所有常見的二維和三維圖表的繪製方法，包括簡單的柱形圖系列、條形圖系列、折線圖系列，地圖系列等。張杰，數據分析-工業設計獅；程序猿；科研狗；簡稱：獅猿狗。
  • Python數據可視化-seaborn Iris鳶尾花數據
    p=5261首先介紹一下Iris鳶尾花數據集，內容摘自百度百科：Iris數據集是常用的分類實驗數據集，由Fisher, 1936收集整理。「Iris也稱鳶尾花卉數據集，是一類多重變量分析的數據集。數據集包含150個數據集，分為3類，每類50個數據，每個數據包含4個屬性。
  • Python數據可視化教程之基礎篇
    經過學習之後，我總結了利用python實現可視化的三個步驟：確定問題，選擇圖形轉換數據，應用函數參數設置，一目了然python中最基本的作圖庫就是matplotlib，是一個最基礎的Python可視化庫，一般都是從matplotlib上手Python數據可視化，
  • Python的數據可視化:對比7種工具包
    Python 的scientific stack（一個介紹Python科學計算包的網站）已經完全成熟，並且有各種各樣用例的庫，包括機器學習（連結：machine learning），數據分析（連結：data analysis）。數據可視化是探索數據和清晰的解釋結果很重要的一部分，但是Python在過去卻相對於其他工具比如R落後一點。
  • Python數據可視化—Seaborn
    今天為大家分享的小技巧是python的可視化畫圖庫Seaborn。相信很多小俠客用過matplotlib、pyecharts等可視化庫，可是為什麼還要介紹Seaborn呢？因為它修復了上述兩個庫的一些缺點，比如Seaborn提供了大量的高級接口和自定義主題，而matplotlib沒有這些接口使得很難確定哪些設置來自定義圖表。
  • 15種獨一無二的創新數據可視化方式
    utm_source=blog&utm_medium=11-data-visualizations-python-r-tableau-d3js2. 用R語言創建BBC式的可視化嚴格意義上來講，這不是可視化，不過也基本上如此。上圖展示了BBC數據團隊發布的可視化融合。
  • Python視頻教程網課編程零基礎入門數據分析網絡爬蟲全套Python...
    1-1課程導學 1-2 數據分析概述 2-1數據倉庫 2-2監視與抓取 2-3填寫,埋點,日誌,計算 2-4 數據學習網站 3-01數據案例介紹 3-02集中趨勢,離中趨勢 3-03數據分布-偏態與峰度 3-04抽樣理論 3-05編碼實現
  • ​數據可視化 | 6個基本可視化Python庫
    數據可視化是使用可視元素（例如圖表，圖形等）表示數據的過程，有助於從數據中獲取有意義的見解。它旨在揭示數據背後的信息，並進一步幫助查看者查看數據中的結構。數據可視化將使任何人只要對數據科學的了解最少，即可獲得科學發現，並幫助人們輕鬆地交流信息。畢竟，一張圖片勝過千言萬語。在本文中我們將介紹一些Python支持的最出色的可視化庫。
  • 使用biopython可視化染色體和基因元件
    在biopython中，通過BiolGraphics子模塊可以對基因組結構進行可視化，支持線性和圈圖兩種可視化方式。其中，基因組結構信息存儲在genebank格式的文件中，首先通過Bio.SeqIO讀取結構信息，然後通過Bio.Graphics模塊進行可視化。
  • 大數據時代,如何培養數據分析思維?|建模|excel|數據分析|python|...
    大數據時代來臨，「數據」熱度飆升，衍生出的行業也受到追捧。據悉，中國大數據行業人才需求2020年將達210萬，未來5年需求量在2000萬人左右。可見，國內數據分析崗錢途與前景並存，想入門該朝陽產業，需要先明白數據分析究竟是什麼？
  • Python爬50W數據,BI做出可視化後,我有了答案
    一次完整的python數據分析流程是怎麼樣的？使用python從網站抓取數據，並將這些數據保存到SQLite資料庫中，然後對數據進行清洗，最後對數據進行數據可視化分析。可是熟悉的人應該知道，python爬取簡單，但是分析起來是很困難的，SQL語句、Pandas和Matplotlib這些十分繁瑣，一般人也不會。
  • 數據統計可視化——python生成詞雲
    有的人會認為它是文本挖掘的可視化、有的人 會給一些用戶打上標籤、有的人則只是想玩玩這酷炫的詞彙……但詞雲產生的視覺效果才是讓人更著迷的地方。詞雲主要是對文本中出現頻率較高的「關鍵詞」進行視覺上的突出，讓人更直觀地看出文本的重點。本文針對上一篇文章
  • Python的可視化工具概述
    >ggplotBokehpygalPlotly在例子用,我將使用pandas操作數據,並啟動其可視化.在大多數情況下使用這些工具不需要pandas,但是我覺得pandas+可視化工具如此普遍,這是最好的起點。
  • 超適合Python小白的乾貨,Python數據分析50個實戰項目
    2.大話NBA | 用數據帶你回顧喬丹的職業生涯3.分析幾十萬條知乎數據，我挖掘出了這些秘密4.用（大）數據全方位解讀電視劇《大秦帝國之崛起》5.以虎嗅網4W+文章的文本挖掘為例，展現數據分析的一整套流程6.從大數據輿情傳播角度看《三生三世十裡桃花》7.北京二手房房價分析8.kaggle數據分析實踐項目練習
  • 量化交易之《Python數據分析》
    Numpy：提供數值數組和函數Scipy：科學計算庫Matplotlib：數據可視化 Numpy 數組的優勢1 Numpy數組能夠運用向量化運算來處理整個數組，而完成同樣的任務，python的列表則通常必須藉助循環語句遍歷列表，並對逐個元素進行相應的處理2 Numpy使用了優化過的C
  • Python 可視化神器--Plotly
    學習Python是做數分析的最基礎的一步，數據分析離不開數據可視化。Python第三方庫中我們最常用的可視化庫是 pandas，matplotlib，pyecharts，當然還有 Tableau，另外最近在學習過程中發現另一款可視化神器-Plotly，它是一款用來做數據分析和可視化的在線平臺，功能非常強大,可以在線繪製很多圖形比如條形圖、散點圖、餅圖、直方圖等等。