背景说明：帕累托系数是用来衡量城市规模位序分布的统计量，回归的目的是计算全国31个省份每个省份每一年的帕累托系数，每一个帕累托系数由某个省份某年中所有城市的观测值回归得到
1.回归方程
ln(rank−1/2)=A−ξ×lnsize+μ
变量说明：size是指城市的人口数据，rank是值城市的位序（desc），A是截距，μ是误差
2.数据说明
输入数据：省份-年份-城市-人口 面板数据（excel文件形式）
期望输出数据：省份-年份-帕累托系数 面板数据（excel文件形式）
3.主要使用的第三方库
statsmodels 是一个 Python 库，用于拟合统计模型、进行统计测试以及统计数据探索和可视化。它是一个非常强大的工具，可以用于众多的统计数据分析任务，从简单的线性回归到复杂的时间序列分析都能胜任
ols回归
import statsmodels.api as sm
# 进行OLS回归
X = sm.add_constant(x)
model = sm.OLS(Y, X) # 生成模型
results = model.fit() # 进行回归并得到回归结果其中x是自变量，add_constant是添加截距的，这是因为常用的统计库都默认没有截距项（只接受x，y两个参数），所以添加截距项的时候是要给自变量加的，然后生成一个新的自变量，本质是给自变量的旁边再加一列截距列
回归结果的获取
results = model.fit()
# 输出模型的系数result.params有结果的各种系数（包括截距和ln_size的系数），
print("截距:", results.params['const']) # 模型的截距
print("ln_size的系数:", results.params['ln_size']) # ln_size的系数
# 输出R方值（模型拟合度的衡量指标）
print("R方值:", results.rsquared)
# 输出调整后的R方值（考虑了自变量数量的拟合度衡量）
print("调整后的R方值:", results.rsquared_adj)
# 输出系数的p值（用于检验系数的统计显著性）
print("截距的p值:", results.pvalues['const'])
print("ln_size的p值:", results.pvalues['ln_size'])
# 输出系数的置信区间
print("ln_size的置信区间:", results.conf_int().loc['ln_size'])
# 输出模型的整体统计摘要
print(results.summary())
# 输出模型残差
print(results.resid)
# 使用模型进行预测（假设有新的数据 new_x）
new_x = [/* 新数据 */]
new_data = sm.add_constant(new_x)
predictions = results.predict(new_data)
print("预测值:", predictions)4.实现代码
import pandas as pd
import numpy as np
import statsmodels.api as sm
def calculate_pareto_coefficient(df):
"""
计算帕累托系数
:param df: DataFrame，包含某个省份某年的所有城市及其人口数据
:return: 帕累托系数
"""
# 按人口数量降序排序并计算排名
df = df.sort_values(by='城市区域人口数量', ascending=False)
df['rank'] = range(1, len(df) + 1)
# 应用帕累托回归公式
df['ln_rank'] = np.log(df['rank'] - 0.5)
df['ln_size'] = np.log(df['城市区域人口数量'])
# 进行OLS回归
X = sm.add_constant(df['ln_size'])
model = sm.OLS(df['ln_rank'], X)
results = model.fit()
# 提取帕累托系数
return -results.params['ln_size']
if __name__ == '__main__':
# 加载数据
data = pd.read_excel('省份_年份_城市_人口_panel_data.xlsx')
# 计算每个省份每年的帕累托系数
pareto_coefficients = data.groupby(['province', '年份']).apply(calculate_pareto_coefficient)
# 将Series转换为DataFrame
pareto_coefficients_df = pareto_coefficients.reset_index()
pareto_coefficients_df.columns = ['Province', 'Year', 'Pareto_Coefficient']
# 将结果保存到新的Excel文件
pareto_coefficients_df.to_excel('pareto_coefficients.xlsx', index=False)
5.代码重点解析
对于代码中的calculate_pareto_coefficient中的回归部分已经解释的比较清楚了，调用函数谁不会啊，最麻烦的还是回归前的数据整理部分。下面详细解释
由于数据是面板数据，所以没有办法直接用dataframe的index和column来直接索引出想要的数据，每一次的回归需要的都是某一省份某一年的所有城市数据，那么我们每一次回归都要筛选出这种维度的数据来，所以要对原始的df进行两次分组，先按省份再按年份，之后每一组都是某个省份某个年份的数据，且分组后每组中包含了匹配的所有数据，形成一个新的df，并不会改变原数据的结构，每组数据都会作为实参传入calculate_pareto_coefficient函数中
# 计算每个省份每年的帕累托系数
pareto_coefficients = data.groupby(['province', '年份']).apply(calculate_pareto_coefficient)
