一直以来，我们都会有一个直观的想法：休息日要比工作日更快乐些，因为不上班嘛，比如下面互联网界流行的一周和一天表情分布图：
再来看看数据的字段和数量，以及是否存在空值。
密集恐惧症的小伙伴估计要被吓尿了！上图根本看不出啥名堂来!
那我们转换一下思路，试试其它的可视化呈现方式~
效果好多了！
上图似乎很好地印证了“人们在深夜不太快乐”的假设，因为图示中可以看出：
从凌晨1点到5点的情绪比率远低于一天中的其他时间点。
让我们再看看一周中的情绪比率分布~
周六和周日的编码是5、6，也就是在图示横轴的最右边，所以上图似乎证明了“人们在周末更快乐”的假设，至少在周六是这样。除非…
我们可以通过根据与上面的2种情景来生成合成数据来检验上面的假设。 先设置生成样本的平均数和样本数。
看起来，早上的情绪更消极些！
然而，如果我们按一周中的每一天(by day-of-the-week)进行细分，该效果就会消失，上面所说的辛普森悖论就出现了。。。
再次强调的是，在构建的数据中，从午夜0点到早上6点的周末没有数据的。
开篇提及的两位论文作者G&M确实在一周中的每一天(day of the week)把数据分开了。他们显示，每小时的趋势是成立的，但有一个需要注意的地方：
在周末，情绪积极的趋势被推迟了两个小时。
更多干货内容，请关注公众号鄙喵的Social Listening与文本挖掘
- Python数据科学实操，包括dataframes、数据可视化和统计分析
- Simpson's Paradox， 其中出现在聚合数据中的结果不会出现在细分数据中
- 多元回归，可确定多个预测因子对结果变量的影响。
1 数据加载和预处理
我们将首先加载数据，它以csv的形式出现。这些数据可以通过爬虫的方式得到，具体操作请问谷歌~ 为了以有用的形式访问这些数据，我们将把它作为dataframe加载。你可以将dataframe视为可以通过编程控制的电子表格(数据的行和列)。 在python中，处理dataframe的标准库是pandas
。让我们使用
pandas
来加载数据。
import pandas as pdimport osdf = pd.read_csv('/home/kesci/input/xps7589/moment_text.txt',delimiter='\t')df.head(3)df.info()
RangeIndex: 99999 entries, 0 to 99998
Data columns (total 4 columns):
hour 99999 non-null int64
pa 99999 non-null float64
na 99999 non-null float64
uid 99999 non-null int64
dtypes: float64(2), int64(2)
memory usage: 3.1 MB需要说明的是，每一行是一条记录，对应一条消息，包含以下字段：
- hour，按每周计算，则最大值为 24*7 - 1 =167 24∗7−1=167" role="presentation" style="box-sizing: border-box; display: inline-block; line-height: normal; word-spacing: normal; overflow-wrap: normal; float: none; direction: ltr; max-width: none; max-height: none; min-width: 0px; min-height: 0px; border-width: 0px; border-style: initial; border-color: initial;">
- pa：正面词汇占比
- na：负面词汇占比
- uid：用户唯一标识符
print("max values:",df.max())print()print("min values:",df.min())max values: hour 167.0
pa 1.0
na 1.0
uid 1956.0
dtype: float64
min values: hour 0.0
pa 0.0
na 0.0
uid 1.0
dtype: float64df['ratio'] = df['pa'] / (df['pa'] + df['na'] + 1e-10)print(df.head(3)) hour pa na uid ratio
0 144 0.083333 0.023810 1 0.777778
1 23 0.090278 0.027778 1 0.764706
2 127 0.200000 0.000000 1 1.0000002 数据描述性可视化
我们要检验以下假设：- 人们周末更快乐
- 人们在深夜不太高兴
作为第一步，让我们先进行浅一些的描述性可视化，也就是各种折线图、散点图之类的。为此，我们将使用matplotlib和seaborn库。
import matplotlib.pyplot as pltimport seaborn as sns%matplotlib inlineplt.scatter(df['hour'],df['ratio'])plt.xlabel('hour')plt.ylabel('sentiment ratio')那就用seaborn绘制一下折线图，变现连续性的变化规律。
sns.lineplot(x='hour',y='ratio',data=df)seaborn
的
lineplot
显示了一周中每个时点的平均情绪。
在默认情况下，误差条是单个标准差。
上面的折线图表明，这些数据之中存在着某种结构，但需要我们细细挖掘出来~
让我们“乘胜追击”，试着梳理出一天中的小时和一周中的天的可视化效果。我们将通过向dataframe添加额外字段来完成此操作。
先对时间，也就是小时时点做处理，得到
时点(time_of_day，0 ~ 23时)
和
天(day_of_week，周一 ~ 周日)
这两个新字段。
df['time_of_day'] = df['hour'] % 24df['day_of_week'] = (df['hour']/24).apply(lambda x : int(x))紧接着，绘制折线图：
sns.lineplot(x='time_of_day',y='ratio',data=df)sns.lineplot(x='day_of_week',y='ratio',data=df);3 辛普森悖论(Simpson's Paradox)
上面欲言又止的以外情形是指“辛普森悖论”，某百科上的解释是：
当人们尝试探究两种变量(比如新生录取率与性别)是否具有相关性的时候，会分别对之进行分组研究。然而，在分组比较中都占优势的一方，在总评中有时反而是失势的一方。该现象于20世纪初就有人讨论，但一直到1951年，E.H.辛普森在他发表的论文中阐述此一现象后，该现象才算正式被描述解释。后来就以他的名字命名此悖论，即辛普森悖论。
总的说来，辛普森悖论就像是欲比赛100场篮球以总胜率评价好坏，于是有人专找高手挑战20 场而胜1场，另外80场找平手挑战而胜40场，结果胜率41%，另一人则专挑高手挑战80场而胜8场，而剩下20场平手打个全胜，结果胜率为28%，比 41%小很多，但仔细观察挑战对象，后者明显较有实力。
量与质是不等价的，无奈的是量比质来得容易量测，所以人们总是习惯用量来评定好坏，而此数据却不是重要的。
除了质与量的迷思之外，辛普森悖论的另外一个启示是：
如果我们在人生的抉择上选择了一条比较难走的路，就得要 有可能不被赏识的领悟，所以这算是怀才不遇这个成语在统计上的诠释。接下来的分析中，你将会看到这个悖论的魔幻之处~
在本次分析中，请小伙伴们考虑以下2个情景：- 平日里，人们通常都不快乐，而且在任何时候都有同样的机会发帖。
- 在周末，人们一般都很快乐，晚上不发帖
我们可以通过根据与上面的2种情景来生成合成数据来检验上面的假设。 先设置生成样本的平均数和样本数。
weekend_mean_ratio = 0.65 # 周末的平均情绪比率weekday_mean_ratio = 0.6 #工作日的平均情绪比率weekday_N = 80000 # 工作日的生成样本数weekend_N = 80000 # 周末的生成样本数首先，让我们生成情绪比率。为了从这个发行版中抽取指定数量的样本，我们将使用numpy库进行数值计算。我们将使用numpy的beta发行版，它生成的值范围介于0到1之间。
import numpy as np正式生成样本：
scale = 2 # 控制样本生成的方差ceweekday_ratios = np.random.beta(scale * weekday_mean_ratio,scale * (1 - weekday_mean_ratio),weekday_N)weekend_ratios = np.random.beta(scale * weekend_mean_ratio,scale * (1 - weekend_mean_ratio),weekend_N)print(f'weekday: {weekday_ratios.mean():0.3f} +- {weekday_ratios.std():0.3f}')print(f'weekend: {weekend_ratios.mean():0.3f} +- {weekend_ratios.std():0.3f}')weekday: 0.600 +- 0.283weekend: 0.650 +- 0.276hours_weekday = np.random.randint(0,24,weekday_N)hours_weekend = np.random.randint(6,24,weekend_N)days_weekday = np.random.randint(0,5,weekday_N)days_weekend = np.random.randint(5,7,weekday_N)df_sim = pd.DataFrame({'time_of_day':np.hstack([hours_weekday,hours_weekend]),'day_of_week':np.hstack([days_weekday,days_weekend]),'ratio':np.hstack([weekday_ratios,weekend_ratios]) } )sns.lineplot(x='time_of_day',y='ratio',data=df_sim)# 添加一列，用来指示星这一天是否为工作日/节假日df_sim['weekday'] = df_sim['day_of_week'].apply(lambda x : x < 5)lineplot = sns.lineplot( x='time_of_day', y='ratio', hue='weekday',data=df_sim )4 多元回归(Multiple regression)
区分这两个潜在的预测因素的一种方法是通过多元回归来验证，即一组 预测目标 (自变量)和 结果 (因变量)之间关系的统计模型。线性回归可以通过以下符号来表示： Y ∼ X1+X2+…+XN. 其中，Y 是因变量，我们将对自变量 X1、X2、…,XN 进行回归。具体地说，我们将估计每个 xi 的系数。在本次分析中，自变量包括是否为工作日(weekday)和一天中的时点(C(time_of_day)[T.n], n表示0 ~ 23时之间的时点)。下面，使用python中的statsmodels做多元回归分析，该包使用与R的glmnet有着类似的语法。
import statsmodels.api as smimport statsmodels.formula.api as smfC(time_of_day)
，指定每个时点都有一个特征的回归。
该公式还包括另一个指示，用于指示一天是否为工作日/休息日。
因变量就是情绪比率。
formula = 'ratio ~ C(time_of_day) + weekday'model = smf.ols(formula='ratio ~ C(time_of_day) + weekday',data=df_sim)results_sim = model.fit()print(results_sim.summary())OLS Regression Results
==============================================================================
Dep. Variable: ratio R-squared: 0.008
Model: OLS Adj. R-squared: 0.008
Method: Least Squares F-statistic: 52.73
Date: Fri, 25 Oct 2019 Prob (F-statistic): 2.74e-251
Time: 07:52:48 Log-Likelihood: -23080.
No. Observations: 160000 AIC: 4.621e+04
Df Residuals: 159975 BIC: 4.646e+04
Df Model: 24
Covariance Type: nonrobust
========================================================================================
coef std err t P>|t| [0.025 0.975]
----------------------------------------------------------------------------------------
Intercept 0.6505 0.005 127.080 0.000 0.640 0.661
C(time_of_day)[T.1] -0.0091 0.007 -1.336 0.182 -0.023 0.004
C(time_of_day)[T.2] -0.0022 0.007 -0.319 0.750 -0.016 0.011
C(time_of_day)[T.3] 0.0062 0.007 0.905 0.366 -0.007 0.020
C(time_of_day)[T.4] -0.0071 0.007 -1.032 0.302 -0.021 0.006
C(time_of_day)[T.5] 0.0084 0.007 1.231 0.218 -0.005 0.022
C(time_of_day)[T.6] 0.0005 0.006 0.078 0.938 -0.011 0.012
C(time_of_day)[T.7] 0.0027 0.006 0.463 0.644 -0.009 0.014
C(time_of_day)[T.8] 0.0002 0.006 0.033 0.974 -0.011 0.012
C(time_of_day)[T.9] -0.0038 0.006 -0.649 0.516 -0.015 0.008
C(time_of_day)[T.10] -0.0019 0.006 -0.316 0.752 -0.013 0.010
C(time_of_day)[T.11] -0.0010 0.006 -0.161 0.872 -0.013 0.011
C(time_of_day)[T.12] -0.0007 0.006 -0.118 0.906 -0.012 0.011
C(time_of_day)[T.13] 0.0017 0.006 0.298 0.766 -0.010 0.013
C(time_of_day)[T.14] 0.0018 0.006 0.312 0.755 -0.010 0.013
C(time_of_day)[T.15] -0.0048 0.006 -0.820 0.412 -0.016 0.007
C(time_of_day)[T.16] 0.0002 0.006 0.033 0.974 -0.011 0.012
C(time_of_day)[T.17] -0.0045 0.006 -0.762 0.446 -0.016 0.007
C(time_of_day)[T.18] -0.0029 0.006 -0.488 0.626 -0.014 0.009
C(time_of_day)[T.19] 0.0021 0.006 0.359 0.720 -0.009 0.014
C(time_of_day)[T.20] -0.0033 0.006 -0.566 0.571 -0.015 0.008
C(time_of_day)[T.21] -0.0064 0.006 -1.079 0.281 -0.018 0.005
C(time_of_day)[T.22] 0.0007 0.006 0.123 0.902 -0.011 0.012
C(time_of_day)[T.23] 0.0029 0.006 0.485 0.628 -0.009 0.014
weekday[T.True] -0.0493 0.002 -32.636 0.000 -0.052 -0.046
==============================================================================
Omnibus: 34097.612 Durbin-Watson: 2.004
Prob(Omnibus): 0.000 Jarque-Bera (JB): 11436.880
Skew: -0.448 Prob(JB): 0.00
Kurtosis: 2.044 Cond. No. 40.7
==============================================================================Warnings:
[1] Standard Errors assume that the covariance matrix of the errors is correctly specified.理解回归结果
a.正如我们所希望的，time_of_day的系数几乎都不显著 ，即p 值大于 0.05。
b.摘要还包括一些关于拟合优度的信息，包括R2、F、AIC和对数似然。这些统计数据有助于决定是否有理由使用更复杂的模型。
回归与统计学意义
time-of-day
变量的 p 值在 p<0.05 的修正值后不会继续存在，但
weekday
变量的 p-值会继续存在。
多重共线性(Multi-collinearity)
C(time_of_day)[T.0]
。
a.这个预测值完全由其他的
time_of_day
变量决定。对于
weekday[T.False]
”和
weekday[T.True]
也是如此。
b.包括这些预测器会使估计问题变得不确定：我们可以在所有
time_of_day预测结果
上加一个常数，然后从
intercept(截距)
项中减去相同的常数。
c.在回归中，这称为多重共线性(multi-collinearity)。当预测因子密切相关时，不可能解释回归系数。
d.由于在本次不是预测型分析，只是对自变量和因变量之间的关系进行解释性分析，这就是为什么密切相关的变量可以毫无顾虑地用于该回归模型的原因。
回归模型可视化
hour_coeffs = [0] # coefficientshour_ses = [0] # standard errorshours = range(1,24)for hour in hours: hour_coeffs.append(results_sim.params[f'C(time_of_day)[T.{hour}]']) hour_ses.append(results_sim.bse[f'C(time_of_day)[T.{hour}]'])对合成数据所反映的规律进行可视化展示：
plt.errorbar([0] + list(hours), hour_coeffs, yerr=hour_ses)plt.hlines(0,xmin=-1,xmax=hours[-1]+1,color='r',alpha=0.5)plt.axis([-1,hours[-1],1.2 * min(hour_coeffs) - max(hour_ses), 1.2 * max(hour_coeffs) + max(hour_ses)])plt.xlabel('hour')plt.ylabel('sentiment coefficient')plt.title('regression results on synthetic data');5 回到真实的数据
最后，让我们对原始数据进行同样的回归分析。df['weekday'] = df['day_of_week'].apply(lambda x : x < 5)model = smf.ols(formula='ratio ~ C(time_of_day) + weekday',data=df)results_real = model.fit()print(results_real.summary())OLS Regression Results
==============================================================================
Dep. Variable: ratio R-squared: 0.002
Model: OLS Adj. R-squared: 0.002
Method: Least Squares F-statistic: 7.621
Date: Fri, 25 Oct 2019 Prob (F-statistic): 2.18e-26
Time: 07:52:58 Log-Likelihood: -49735.
No. Observations: 99999 AIC: 9.952e+04
Df Residuals: 99974 BIC: 9.976e+04
Df Model: 24
Covariance Type: nonrobust
========================================================================================
coef std err t P>|t| [0.025 0.975]
----------------------------------------------------------------------------------------
Intercept 0.6077 0.007 87.974 0.000 0.594 0.621
C(time_of_day)[T.1] -0.0291 0.010 -2.853 0.004 -0.049 -0.009
C(time_of_day)[T.2] -0.0340 0.011 -3.075 0.002 -0.056 -0.012
C(time_of_day)[T.3] -0.0423 0.012 -3.449 0.001 -0.066 -0.018
C(time_of_day)[T.4] -0.0303 0.013 -2.405 0.016 -0.055 -0.006
C(time_of_day)[T.5] -0.0346 0.013 -2.740 0.006 -0.059 -0.010
C(time_of_day)[T.6] -0.0034 0.012 -0.290 0.772 -0.026 0.019
C(time_of_day)[T.7] 0.0067 0.010 0.663 0.507 -0.013 0.027
C(time_of_day)[T.8] 0.0197 0.009 2.120 0.034 0.001 0.038
C(time_of_day)[T.9] 0.0209 0.009 2.367 0.018 0.004 0.038
C(time_of_day)[T.10] 0.0298 0.009 3.430 0.001 0.013 0.047
C(time_of_day)[T.11] 0.0353 0.009 4.097 0.000 0.018 0.052
C(time_of_day)[T.12] 0.0206 0.009 2.390 0.017 0.004 0.037
C(time_of_day)[T.13] 0.0183 0.009 2.131 0.033 0.001 0.035
C(time_of_day)[T.14] 0.0168 0.009 1.958 0.050 -2e-05 0.034
C(time_of_day)[T.15] 0.0140 0.008 1.653 0.098 -0.003 0.031
C(time_of_day)[T.16] 0.0048 0.008 0.568 0.570 -0.012 0.021
C(time_of_day)[T.17] 0.0098 0.009 1.155 0.248 -0.007 0.027
C(time_of_day)[T.18] -0.0031 0.009 -0.366 0.714 -0.020 0.014
C(time_of_day)[T.19] -0.0036 0.009 -0.418 0.676 -0.020 0.013
C(time_of_day)[T.20] 0.0116 0.009 1.363 0.173 -0.005 0.028
C(time_of_day)[T.21] 0.0200 0.008 2.354 0.019 0.003 0.037
C(time_of_day)[T.22] 0.0083 0.009 0.974 0.330 -0.008 0.025
C(time_of_day)[T.23] 0.0020 0.009 0.225 0.822 -0.015 0.019
weekday[T.True] -0.0075 0.003 -2.682 0.007 -0.013 -0.002
==============================================================================
Omnibus: 298895.639 Durbin-Watson: 1.857
Prob(Omnibus): 0.000 Jarque-Bera (JB): 11373.330
Skew: -0.517 Prob(JB): 0.00
Kurtosis: 1.711 Cond. No. 33.8
==============================================================================
Warnings:
[1] Standard Errors assume that the covariance matrix of the errors is correctly specified.现在，几个时点的系数有 显著的p值，即使是在Bonferroni校正之后，这些p值也是显著的。这表明这几个时点的情绪比率存在特殊性，也就是规律性，如1 ~ 5点的情绪较负面。至于影响的力度大小，可以看系数，也就是上面结果中的coef 。
有点意思！
hours = range(1,24)hour_coeffs = [results_real.params[f'C(time_of_day)[T.{hour}]'] for hour in hours]hour_ses = [results_real.bse[f'C(time_of_day)[T.{hour}]'] for hour in hours]plt.errorbar(list(hours), hour_coeffs, yerr=hour_ses)plt.hlines(0,xmin=-1,xmax=hours[-1]+1,color='r',alpha=0.5)plt.axis([-1,hours[-1],1.2 * min(hour_coeffs) - max(hour_ses), 1.2 * max(hour_coeffs) + max(hour_ses)])plt.xlabel('hour')plt.ylabel('sentiment coefficient')plt.title('Regression Results on Real Data');可以做的一些拓展分析
- 描述另一种辛普森悖论(Simpson's paradox)，即每个发布消息的用户不受星期几或一天中的时间影响，但仍然存在总体水平(population-level)的影响。
- 在提供的数据中使用多元回归来控制此问题。在statsmodel的回归中添加C(uid)可能会抛出MemoryError，因此您必须考虑其他方法来进行回归。
- 你可能会发现
scipy.sparse.csr_matrix很有用。
更多干货内容，请关注公众号鄙喵的Social Listening与文本挖掘。
451
被折叠的 条评论
为什么被折叠?
到【灌水乐园】发言
