伦敦犯罪数据分析
我在研究伦敦的犯罪数据集。分析和可视化比我预期的要好得多，所以我想我还不如写一篇关于它的文章。我在试验 bigquery 时获得了这些数据，但是你也可以从 Kaggle 上的这里下载。让我们从加载库和数据开始。
library(tidyverse) ## For data wrangling and visualization
library(lubridate) ## To work with dates
library(ggpubr) ## Extra visualizations and themes
library(patchwork) ## Patch visualizations together
library(hrbrthemes)## extra themes and formatting
library(scales) ## For formatting numeric variables
library(tidytext) ## Reordering within facets in ggplot2
library(pier) ## Make interactive piecharts in R
library(ggalt) ## Extra visualizationslondon_crimes <- read_csv("C:\\Users\\ACER\\Downloads\\bq-results-20190830-212006-dh1ldd1zjzfc.csv",trim_ws = TRUE) %>%
mutate(Date_Column = dmy(Date_Column))
加载完库之后，我使用了 read_csv 函数，并给了它一个文件在我电脑上的路径，设置 trim_ws 为 TRUE 以防有任何空白并删除它，然后使用 mutate 函数将 Date_Column 更改为 Date 类型，使用的是Lubridate包中的 dmy 函数。
为了得到一个概览，我将想象一下这些年来整体犯罪事件是如何变化的。让我们浏览一下代码:
london_crimes %>%
group_by(Year = floor_date(Date_Column,unit = "year")) %>%
summarise(Incidents=sum(total_incidents,na.rm = TRUE)) %>%
ungroup() %>% mutate(pct_change= (Incidents-lag(Incidents))/lag(Incidents),
pct_change=replace_na(pct_change,0)) %>%
ggplot(aes(Year,Incidents))+
geom_bar(stat="identity",fill="firebrick",color="black")+
geom_line(color="steelblue",size=1.5,linetype="dashed")+
geom_text(aes(label=percent(pct_change)),vjust=-1,color="black",face="bold")+ geom_text(aes(label=comma(Incidents)),vjust=1,fontface="bold",color="white")+
scale_y_comma(expand = c(0,0),limits = c(0,800000))+
scale_x_date(breaks = "year",date_labels ="%Y")+
theme_classic()+
labs(title = "Total Incidents over the years")
- 为了获得年度分析，我们按日期列分组，并使用 floor_date 函数将其四舍五入，然后使用summary函数获得每年的事件总数
- 每次 group_by 操作后，取消分组总是一个好的做法，以删除分组上下文
- 然后，我们使用 mutate 函数创建一个新列— pct_change ，它基本上是使用 lag 函数计算每年之间的百分比变化
- 整理好数据后，我们使用 ggplot2 制作一个标准的条形图。由于 ggplot2 不允许组合图表显示百分比增长线图，我只是使用了 geom_line 函数来绘制趋势线，并用百分比变化值对其进行注释
如果我们只是比较 2008 年和 2016 年，你会认为犯罪事件几乎没有变化，但实际上这是 2014 年的最低水平，此后一直在增长。现在让我们在区一级分析我们的数据。
让我们来看一个简单的柱状图，看看区一级的事故总数。如果你熟悉 Tidyverse 元包，那么下面这些应该不难理解。
london_crimes %>%
group_by(borough) %>%
summarise(Incidents=sum(total_incidents)) %>%
ungroup() %>%
ggplot(aes(reorder(borough,Incidents),Incidents))+
geom_bar(stat = "identity",aes(fill=borough),color="black")+
coord_flip()+
scale_y_comma()+
geom_text(aes(label=comma(Incidents)),hjust=1)+
theme_classic()+
theme(legend.position = "none")+
labs(x=" ",y=" ",title = "Total Incidents for boroughs from 2008-2016 ")
我不熟悉伦敦的地理，但 Westminister 遥遥领先。我谷歌了一下，发现 Westminister 是伦敦人口最密集的地区之一。
这些年来这些区的犯罪率有什么变化？
为此，我将 2016 年的数据与 2014 年的数据进行比较，而不是与 2008 年进行比较，因为 2014 年的犯罪事件最低，与 2016 年进行比较可能会显示一些有趣的变化。
london_crimes %>%
group_by(borough) %>% summarise(Incidents_2014=sum(total_incidents[year(Date_Column)==2014]), Incidents_2016=sum(total_incidents[year(Date_Column)==2016])) %>%
ungroup() %>%
mutate(Pct_Change=(Incidents_2016-Incidents_2014)/Incidents_2016) %>%
ggplot(aes(reorder(borough,Pct_Change),Pct_Change))+
geom_bar(stat = "identity",aes(fill=borough),color="black")+
coord_flip()+
scale_y_continuous(labels = percent_format())+
geom_text(aes(label=percent(Pct_Change)),hjust=1)+
theme_classic()+
theme(legend.position = "none")+
labs(x=" ",y="Percentage Change ",title = "Percentage Change in Incidents from 2014-2016")
- 我们再次按区对数据进行分组，并计算两个条件和——Incidents _ 2004 和 Incidents_2016，并计算两者之间的百分比差异，以查看每个区的变化
- 然后，像以前一样，按照标准程序绘制柱状图，根据从最高到最低的百分比变化重新排序
- 伦敦金融城的最高比例有点失真，因为它的犯罪事件数量很少
- 令人惊讶的是，Westminister 是倒数第二名，尽管它在所有事故中排名第一。至少，它没有经历相当大的增长。
我想知道每一个行政区在哪一年经历了有史以来最高的犯罪事件数量，以及哪一年有最多的行政区经历了最高的犯罪事件数量。
london_crimes %>%
group_by(borough,Yearly=floor_date(Date_Column,unit = "year")) %>%
summarise(Incidents=sum(total_incidents)) %>%
ungroup() %>%
group_by(borough) %>%
filter(Incidents==max(Incidents)) %>%
ungroup() %>%
mutate(borough_max_year = paste0(borough,"-","(",year(Yearly),")")) %>%
ggplot(aes(reorder(borough_max_year,Incidents),Incidents))+
geom_col(aes(fill=borough_max_year),color="black")+
scale_y_comma()+
coord_flip()+
geom_text(aes(label=comma(Incidents)),hjust=1)+
theme_classic()+
theme(legend.position = "none")+
labs(title = "Max Incidents for each Borough",x="Borough and year of max Incidents ",y="Incidents")
- 我们首先按照年份和行政区进行分组，然后根据行政区和相应年份的每个组合获得事故总数
- 因为我想根据事件筛选每个区的最大年份，所以我按区对数据进行分组，然后筛选每个区的事件值等于其最大值的地方。这也将显示该值的年份
- 为了直观起见，我创建了一个新列— borough_max_year ，它基本上是使用 paste0 函数将year和 Incidents 列串联起来的版本。剩下的是我用来可视化柱状图的相同代码。
现在，我们可以看到每一个区在哪一年经历了最高数量的犯罪事件，但我现在将制作另一个可视化效果，以查看哪一年最多的区经历了最高数量的犯罪事件。
london_crimes %>%
group_by(borough,Yearly=year(floor_date(Date_Column,unit = "year"))) %>%
summarise(Incidents=sum(total_incidents)) %>%
ungroup() %>%
group_by(borough) %>%
filter(Incidents==max(Incidents)) %>%
ungroup() %>%
count(Yearly,sort = TRUE,name = "Boroughs_with_max_incidents") %>% ggplot(aes(reorder(Yearly,Boroughs_with_max_incidents),Boroughs_with_max_incidents))+
ggalt::geom_lollipop(aes(color=factor(Yearly)),point.size = 10)+
scale_y_comma()+
coord_flip()+
geom_text(aes(label=Boroughs_with_max_incidents),hjust=0.5,color="white",fontface="bold")+
theme_classic()+
theme(legend.position = "none")+
labs(title = "Number of boroughs that experienced their highest ever incidents",x="Year",y=" ")
- 使用和以前一样的代码，我增加了一个额外的步骤，计算每年的行政区数量
- 厌倦了条形图，我使用了 ggalt 包中的 geom_lollipop 函数来绘制棒棒糖图
- 因此，2016 年拥有最多的行政区，尽管 2016 年在 2008 年和 2012 年之前的犯罪事件总数方面排名第三
我现在要检查每个区的犯罪事件的总体趋势。
london_crimes %>%
group_by(Year=floor_date(Date_Column,unit = "year"),borough) %>%
summarise(Incidents=sum(total_incidents,na.rm = TRUE)) %>%
ggplot(aes(Year,Incidents))+
geom_line(aes(color=borough),size=0.75,style="--")+
theme_pubclean()+
scale_y_comma()+
expand_limits(y=0)+
facet_wrap(~borough,scales = "free")+
labs(y="Total Incidents",x=" ")+
theme(legend.position = "none",strip.background = element_rect(fill="firebrick"),strip.text=element_text(color = "white",face="bold"))
我原以为这将是一个非常混乱的图表，但结果比我预期的要好很多。从 2014 年到 2016 年，我们可以大致看到大多数行政区的总体增长趋势。
让我们来看看每个主要类别在犯罪事件总数中所占的比例，以及它们的长期趋势。
**# Piechart**london_crimes %>%
group_by(major_category) %>%
summarise(Incidents=sum(total_incidents,na.rm = TRUE)) %>%
ungroup() %>%
mutate(color=RColorBrewer::brewer.pal(9, 'Spectral')) %>%
select(label=major_category,value=Incidents,color) %>%
pier() %>%
pie.size(inner=70, outer=100, width = 600, height = 450) %>%
pie.header(text='Crimes', font='Impact', location='pie-center') %>%
pie.subtitle(text='by major category')%>%
pie.tooltips()**# Trend over time**london_crimes %>%
group_by(Yearly=floor_date(Date_Column,unit = "year"),major_category) %>%
summarise(Incidents=sum(total_incidents,na.rm = TRUE)) %>%
ggplot(aes(Yearly,Incidents))+
geom_line(aes(color=major_category),size=0.75)+
theme_pubclean()+
scale_y_comma()+
expand_limits(y=0)+ # making sure that the yaxis for every facet starts at 0 otherwise the trend may look misrepresentitive
facet_wrap(~major_category,scales = "free")+
labs(y="Total Incidents",x=" ")+
theme(legend.position = "none",strip.background = element_rect(fill="firebrick"),strip.text=element_text(color = "white",face="bold"))
- 我在这里使用了 pier 包来绘制一个交互式饼图，因为我总是发现在 ggplot2 中绘制像样的饼图比实际需要的要复杂一些
- pier 要求栏目只需标注“标签”、“数值”和“颜色”
- 趋势图的代码与行政区的代码相同。唯一的区别是我用 major_category 变量替换了 borough 变量
- 请注意，我使用了 expand_limits(y=0) 参数来确保每个线图中的 y 轴都从 0 开始，否则，它们看起来会有点失真
- 从 2014 年到 2016 年，盗窃和处理以及针对人身的暴力等主要贡献者出现了显著增长
- 似乎奇怪的是，2008 年后就没有发生过性犯罪和欺诈或伪造。可能是因为这些罪行后来被记录在不同的类别下
现在我想看看不同的行政区在每个主要类别中的排名。我将只关注前 40%的地区，因为在一个多面图中可视化所有的行政区看起来非常混乱。
london_crimes %>%
group_by(major_category,borough) %>%
summarise(Incidents=sum(total_incidents)) %>%
filter(Incidents>quantile(Incidents,0.6)) %>% #filtering for only the top 40% to make the plot more readable
ungroup() %>%
mutate(borough=reorder_within(borough,Incidents,major_category)) %>% ggplot(aes(borough,Incidents))+
geom_bar(aes(fill=borough),stat = "identity",color="black")+
coord_flip()+
scale_x_reordered()+
scale_y_comma()+
geom_text(aes(label=comma(Incidents)),hjust=1)+
theme_classic()+
theme(legend.position = "none")+
facet_wrap(~major_category,scales ="free")+
labs(x=" ",y=" ",title = "Top 40% Boroughs by Incidence",subtitle = "Segmented by Major Category")+
theme(legend.position = "none",strip.background = element_rect(fill="firebrick"),strip.text=element_text(color = "white",face="bold"))
- 在分组和总结之后，我保留分组上下文，使用分位数函数过滤每个主要类别下前 40%的行政区
- 然后，我使用 tidytext 包中的 reorder_within 函数，根据犯罪事件对每个主要类别下的行政区进行重新排序。这样做的原因是为了在每个主要类别方面内对条形图进行排序
另外，我想看看每个大类中哪些是主要的小类。
london_crimes %>%
group_by(major_category,minor_category) %>%
summarise(Incidents=sum(total_incidents)) %>%
ungroup() %>%
mutate(minor_category=reorder_within(minor_category,Incidents,major_category)) %>%
ggplot(aes(minor_category,Incidents))+
geom_bar(aes(fill=minor_category),stat = "identity",color="black")+
coord_flip()+
scale_x_reordered()+
scale_y_comma()+
geom_text(aes(label=comma(Incidents)),hjust=1)+
theme_classic()+
theme(legend.position = "none")+
facet_wrap(~major_category,scales ="free")+
labs(x=" ",y=" ",title = "Incidents by Minor Category",subtitle = "Segmented by Major Category")+
theme(legend.position = "none",strip.background = element_rect(fill="firebrick"),strip.text=element_text(color = "white",face="bold"))
- 所以大部分都是骚扰、袭击和盗窃案件。我猜在 2008 年后，性犯罪被归入对人身的暴力。
- 看起来商业和个人财产相对来说被抢劫了很多。虽然盗窃和处理有很多潜在的重叠，但是抢劫和入室盗窃。
结论
尽管数据集性质严峻，但我在探索它的过程中获得了很多乐趣。也有可能进行更多的探索性分析。你可以对伦敦最安全的行政区进行分析，或者深入挖掘小类事故。
分析纽约市的交通方式
这项工作由坎贝尔·韦弗、吴宇舟和雷汉·拉苏尔完成，是康奈尔理工大学**数据科学课程的一部分。**
简介
对大都会运输管理局(MTA)来说，预测给定纽约市地铁站的乘客数量是一个关键挑战。
因此，我们旨在解决以下挑战:给定一周中的某一天、一个站点位置和附加数据(例如天气、特殊事件)，他们预计会有多少人？
这种分析和预测的目的是让 MTA 为不可预见的情况做好准备。例如，在罗斯福岛最近举行的樱花节期间，地铁站出现了大拥堵，这使许多人面临窒息或踩踏的危险。
使用我们的机器学习模型，MTA 将能够预测不同天气条件下每个地铁站的预期人数，并考虑到特殊事件。这将允许他们先发制人的安排，以避免任何危险的情况。
这项工作侧重于 2018 年的数据。
分析
数据选择
****地铁数据来自 MTA 十字转门数据。尽管这些数据比我们需要的更复杂，但是我们能够执行彻底的数据清理，使其对我们的分析有用。
该数据由所有纽约地铁站十字转门的入口和出口计数器组成，每 4 小时采集一次。
我们使用这个数据集来分析天气和特殊事件对地铁使用的影响。
数据集大小:1000 万行，11 列
从谷歌地图 API 获得的位置数据。这被用来提取从 MTA 数据获得的台站的纬度和经度，因为这些数据只包含台站名称，而不包含它们相应的位置坐标。
从城市自行车出行历史中获得的城市自行车数据。数据集由 2018 年期间每次骑自行车的信息组成，包括行程持续时间、开始&结束时间、日期和车站、车站位置(纬度/经度)和一些用户信息。
我们使用这个数据集来分析天气和特殊事件对自行车骑行的影响。
数据集大小:1100 万行，7 列
****从国家环境信息中心获得的天气数据。数据集由 2018 年纽约市不同站点的气候数据组成，包括站点位置(纬度/经度)和各种气候参数，例如，温度、降水量和降雪量。
数据集大小:28000 行，33 列
活动数据:吸引大量人群(从而导致地铁乘客高峰)的特殊活动有许多不同的形式。因此，很难找到包含所有这些事件的历史信息的现有数据集。在我们分析事件对交通数据的影响时，我们集中在 161st St. Yankee Stadium 站的乘客量和 Yankee 主场时间表之间的关系。为了获得扬基队的赛程信息，我们使用开源 python 包装器来访问 MLB 的 API。
数据清理
****地铁数据使用以下步骤进行清理:
- 每天每个十字转门的总进出数据
- 每天每个站点的总进出数据
- 使用 z 得分> 3 移除异常值
这产生了一个干净的数据集，包含 136000 行，比我们开始时的 1000 万行减少了很多。
****位置数据用于使用站名填写每个站的位置坐标(纬度/经度)。
****CitiBike 数据非常大，因为每行仅代表一次骑行记录。为了节省时间和空间，我们把数据过一遍，收集每个站的数据。我们将站名、id 和位置存储到字典中以备后用。此后，我们通过计算每日乘车人数和平均乘车时长来压缩数据。因此，我们获得了每天每个站点的自行车使用数据。
****天气数据使用以下步骤进行清理:
- 仅提取纽约市数据
- 填写 TAVG 栏
- 将缺失值填写为 0
- 小写列保持一致
这产生了一个干净的数据集，从 28000 行减少到 365 行，一年中每天一行。
数据分析
下雨会影响地铁的使用吗？
distribution of subway usage for rainy vs non-rainy days
零假设(H0):
雨天和非雨天的地铁使用量没有区别
交替假设(HA):
雨天和非雨天的地铁使用量是有区别的
每个类别有两个分布，我们猜测这与周末和工作日有关。我们用下面的图来证实。
看起来有两种不同的多雨和不多雨的分布，但是我们不能确定这种差异是偶然的还是有统计学意义的。因此，我们对该数据执行非参数 Mann Whitney U-test ，因为样本是独立且非正态分布的。
p 值出来是 0.0218，小于我们 p 临界值的 0.05。因此，我们可以拒绝零假设，并得出结论雨天和非雨天的地铁使用存在差异。****
我们通过绘制雨天和非雨天的地铁使用热图来证实这种差异。
下雨会影响 CitiBike 的使用吗？这是一个更简单的问题，因为人们本能地不太可能在雨天骑自行车。我们通过绘制雨天和非雨天的 CitiBike 使用热图证实了这一点。
地铁和城市自行车的乘客量每月都有变化吗？为了回答这个问题，我们绘制了 subway 和 CitiBike 每月的使用情况。
看起来在炎热的月份，地铁乘客量是最低的，因为人们在那些日子里避免乘坐地铁，而更喜欢骑自行车或步行。
特殊事件对地铁使用有何影响？在体育比赛或音乐会等受欢迎的活动当天，靠近活动场地的地铁站的乘客量会增加，这是很直观的。在许多情况下，这种影响很难量化，因为许多地铁站可以服务于同一个场地，从而将这种影响分散到不同的交通枢纽。此外，许多一次性事件，如音乐会，只产生一天的数据，并且很少重复，留给我们的样本量很小，无法建立我们的直觉。最后，收集纽约市许多不同类型的活动和场所的历史数据是一项艰巨的任务。我们选择通过关注 161 ST . Yankee Stadium(Yankee Stadium 附近唯一方便的地铁站)Yankee 主场比赛对地铁乘客量的影响来分析乘客量和事件之间的关系。
通过绘制该车站一年的乘客数据，一些事情立即变得显而易见。首先，有一个规律周期，即工作日的客流量较高(每天约 10，000 人)，周末的客流量较低(约 5，000 人)。突然，大约在 4 月 2 日，在洋基主场揭幕战的当天，每天的乘客量开始变成穿插着超过 15，000 名乘客的常规高峰。这种情况会持续到洋基的球季在 10 月结束，届时上座率会回复到先前的水平。
通过把我们的电视台数据分成有洋基主场比赛的日子和没有主场比赛的日子，我们立刻看到了不同。当有主场比赛时，平均每天的观众人数是 18，000，当没有主场比赛时，观众人数下降到 10，000。通过观察当有洋基主场比赛时和没有比赛时车站乘客的分布，很明显这些来自不同的分布。
同样，我们分析了 2018 年罗斯福岛地铁交通。2018 年 4 月 21 日地铁使用量有一个明显的峰值，对应的是罗斯福岛樱花节。
预测模型
我们使用了两种机器学习模型来根据位置和天气信息预测车站的客流量:
- 线性回归
- 随机森林
我们使用以下测量技术来评估性能:
- 平均绝对误差(MAE) —误差绝对值的平均值
- 均方根误差(RMSE) —均方误差的平方根
Performance of different models
从表中我们可以看出，与线性回归模型相比，随机森林模型具有更好的性能。由于站点位置对客流量的影响不能用线性表示，因此线性模型很难准确预测。
为了可视化这些预测，我们创建了预测值与实际值的热图。我们可以清楚地看到，随机森林模型优于线性回归模型:
我们还研究了不同特征对乘客量的重要性。
****
经度和纬度对日客流量影响最大。这是有道理的，因为不同的社区有不同的人口，这是公共交通使用的最重要的因素。
平均温度表明了一个普遍的天气状况和一年中的季节。所以也是一个重要的特点。
在工作日和周末，运输使用情况不同，显示为一周中的某一天。
降水也会影响乘客量，但与其他因素相比，它的影响较小。
进一步的改进
可以采取以下步骤来改进这项工作:
- 使用具有更丰富特征的更高级的机器学习模型，例如，我们可以将事件信息作为特征的一部分，以获得更好的预测结果
- 包括出租车和优步/Lyft 数据，以分析它们如何受到天气和事件的影响
- 使用多年的数据捕捉年度趋势
用 Python 分析电影对话速度
当大卫·芬奇在执导《社交网络》时，他带着秒表去了阿伦·索尔金的家，记录了索尔金阅读剧本的时间。然后，他记下了每个场景的时间，并在拍摄这部杰作时严格遵守这些场景时间。听到这个轶事后，我决定自己从分析的角度来看看这部电影的对话。
为了开始我的调查，我首先收集数据。我发现了一个叫做 YIFY 字幕(https://yts-subs.com/)的网站，它提供电影编辑用来给电影添加字幕的文件。这些特定的文件称为 SRT 文件，甚至有自己独特的。srt 扩展。看一些例子。srt 文件，很明显它们遵循一种严格的格式，包含:
1。字幕编号
2。格式为 START_TIME → END_TIME 的场景时间，时间格式为 hh:mm:ss:ffff
3。副标题的文本
我需要删除文本中的一些元素(html 标签、圆括号和方括号中的场景描述符、音乐符号),所以我编写了函数 clean_text 来帮助清理文件。
def clean_text(line):
clean = line.strip().replace('</i>', '').replace('<i>','').replace('- ', '')
clean = re.sub('\(.*\)','', clean)
clean = re.sub('\[.*\]','', clean)
clean = re.sub('♪', '', clean)
return clean
通过几个辅助函数，我遍历了文件中的每一行，获得了以秒为单位的累积电影时间和说话的字数。解析每个的最终输出。srt 文件是一个具有两列的数据帧，一列用于累计时间，一列用于累计说出的单词。
为了有一个比较《社交网络》的电影基础，我还收集了一些波士顿同行的数据:
善意狩猎
狂热程度
战士
消失的婴儿
爱国者日
逝者
乡野圣徒
神秘河
小镇
所有的。srt 文件被解析，数据在干净、统一的数据帧中，然后我继续进行可视化。随着时间的推移，绘制每部电影的字数突出了社交网络，因为它具有最多的单词和最陡的斜率(对应于每秒的单词数)。
据传，芬奇花了 99 次拍摄才对《社交网络》的开场感到满意。将窗口缩小到前 600 秒会进一步提高对话速度。
利用每部电影的最终长度和字数，我还计算了每秒的平均字数。
最后，我又添加了 15 部电影来增加样本量，并绘制在下面。即使像《超级坏蛋》和《婚礼傲客》这样有大量对话和快嘴演员的电影，社交网络仍然保持领先。
感谢阅读！对于那些感兴趣的人，我已经把我用来进行这种分析的笔记本和。我的 GitHub 页面上的 srt 文件(https://github.com/eonofrey/movie_word_counts)。谢了。要获得所有媒体文章的完整访问权限，请点击此处！
用机器学习分析我的减肥
个人健康
我如何用 Python 从头开始构建一个逻辑回归分类器来预测我的体重减轻
要看我为这个项目写的代码，可以看看它的 Github 回购
背景
我在 2018 年初开始了我的减肥之旅，遵循着人们常说的“减肥=饮食+运动”的建议。在饮食方面，我开始跟踪我每天的食物消耗量(使用食物秤并通过 Loseit 应用程序记录卡路里)。在锻炼方面，我开始遵循“沙发到 5 公里”计划，到目前为止已经完成了四个 5 公里，一个 10K，还有几周前的一个半程马拉松。最后，每天早上，我醒来后马上称体重，并在同一个 Loseit 应用程序中记录我的体重。
问题
任何试图减肥的人都不可避免地会遇到减肥平台期，最初快速的减肥开始放缓。就我个人而言，我在五月初的假期后达到了一个主要的平台期。在那之后，我放弃了追踪我的进度，持续了近三个月。只有在我用一个 Amazfit Bip 取代了快要报废的 Pebble 智能手表后，我才重新获得了一些按下恢复按钮的动力，部分原因是我可以用新手表开始记录我的步数。然而，我的体重持续稳定，在经历了两个月令人沮丧的体重波动后(见下面的虚线区域)，我完全停止了追踪我的体重和卡路里。那是去年的 11 月。
Left: Progress report of my weight in 2018. Right: time window for my data (when steps were tracked)
现在是 2019 年，随着 Tet(越南新年)最近接近尾声，我决定更仔细地研究在这两个月期间收集的数据，希望发现我的减肥和跟踪的卡路里/步数之间的有趣关系，以便我可以制定比以前更有效的减肥计划。
数据
数据收集
Calories (budget & consumed) and steps for September 6th, 2018
**卡路里:**从我的 Loseit 账户导出为 CSV 文件。出于某种奇怪的原因，Loseit 只允许一次导出一周的卡路里数据，但是将它们连接在一起是快速 Python 脚本所不能做到的。每个日期都有我当天记录的所有食物的相应卡路里计数，以及应用程序根据我当天的体重和我最初指定的减肥目标(每周减肥 0.75 公斤)为我计算的卡路里“预算”。
**步骤:**我的 Amazfit Bip 智能手表的 Android 应用程序不允许数据导出，除非有人使用第三方工具使用一些令人毛骨悚然的变通方法。因此，获得我的步骤数据的最快方法是手动滚动手机上的几十个日期(在这两个月内)，并将每个日期的步骤输入到一个 CSV 文件中。一点也不优雅，但是，嘿，这很有效！
**体重:**谢天谢地，Loseit 网站允许我将我所有记录的体重——以及它们的拍摄日期——导出为一个 CSV 文件。
在按日期将这 3 个数据源连接在一起后，我最终得到了这两个月中只有 46 个日期的卡路里+步数+体重数据。显然，我在相当多的日子里(出于显而易见的原因，很多日子是周末)忽略了记录这三个数据中的至少一个。
数据转换
根据这三个原始数据源，计算三个附加数据字段:
剩余=消耗的卡路里-卡路里预算。正剩余意味着我吃了比那天允许的预算更多的卡路里，反之亦然。我选择使用卡路里盈余而不是消耗的原始卡路里，因为应用程序的卡路里预算会随着我体重的上升和下降而自然变化，所以卡路里盈余(考虑了上述预算)比卡路里本身更能准确地衡量我的饮食习惯。
增重=明天体重-今天体重。给定日期的正体重增加意味着我在那一天增加了体重(咄！)，反之亦然。
**体重增加状态:**正的体重增加标记为 1，负的或零体重增加标记为 0。我决定使用二元体重增加状态——无论我是否增加了体重——而不是更具体的体重增加量，因为纠结于增加 0.5 公斤还是 0.3 公斤是非常适得其反的，尤其是因为体重增加量会受到除卡路里和步骤之外的许多因素的影响(如饮水量、进食时间等)。
First few rows of my final data table
数据可视化
当绘制卡路里剩余量和每天的步数与我当天体重增加状况的关系图时(见下面的前两个图)，这两个据说很重要的因素如何预测我是否会增加体重似乎没有明显的模式。
Left & middle: plot of weight gain status against calories surplus and steps separately. Right: plot of weight gain status against surplus and steps together
然而，当绘制卡路里剩余量和步数时(上图右侧)，非常有趣的模式出现了！例如，我马上就能看出，根据我的步数，有两组截然不同的日期:那些低于 5000 步的日期(我的“基线”懒惰日)整齐地围绕着一条水平线，而那些高于 5000 步的日期(我的活跃日)，这主要归功于我的跑步。就卡路里过剩而言，有三个主要观察结果:
1。在我懒惰的日子里，如果我吃的超过了我的卡路里预算限制，第二天当我走上体重秤时，这将是个坏消息。有一些奇迹，比如有一天我吃了超过限量的 1500 卡路里，但第二天体重并没有增加，但这种情况很少。
2。另一方面，如果我在懒惰的日子里吃得低于我的卡路里预算限制，我也不完全清楚:有好几天我是个好孩子，吃我的沙拉(打个比方)，但仍然增加了体重。这表明我应该在懒惰的日子里更加保守地饮食。
3。然而，在我活跃的日子里，我似乎可以吃得比我的卡路里限量还多，因为有几天我吃得比限量还多，但体重却没有增加。
这些观察表明，我应该考虑我的日常活动(以步数的形式)，而不是使用 Loseit 应用程序的默认限制。
- 例如，在上面的剩余步骤图上，我可以画一条直线，将我体重增加的天数(红色)和体重减少的天数(绿色)很大程度上分开。
- 然后，我可以使用该线性边界来告知我的减肥策略，即我如何保持在边界的减肥一侧，而不是另一侧。
- 用机器学习的说法，这相当于建立一个线性分类器对我的数据进行二进制分类(对体重增加和体重减少进行分类)。
使用哪个线性分类器？
从数据构建线性分类器有几种常见的方法，如逻辑回归、[线性判别分析](http://linear discriminant analysis)或支持向量机。但是，对于这个项目，我将使用逻辑回归，因为:
- 这很容易解释:例如，分类边界的方程可以很容易地从逻辑回归系数中获得。
- 易于实现:这是一个非常重要的原因，因为这个项目的另一个目标是让我实现一个机器学习项目，而不使用预先存在的库(如 scikit-learn)。正如我将在后面展示的，逻辑回归学习算法的核心只需要 3 行代码就可以完成！
选择模型后，让我们看看如何使用我的数据实现和解释逻辑回归。然而，这需要我们回顾一些关于逻辑回归分类器如何工作以及如何从数据中学习的数学知识。
数学评论
请注意，这个数学复习主要是为了建立一个通用的符号，以便在从数学方程实现算法之前，我们在同一页上。如果你想了解这些方程是如何推导出来的，以及它们背后的直觉，我在这篇博文的末尾链接了更多的资源，这些资源可以比我做得更好。
根据特征预测体重增加的概率
在逻辑回归中，我在某一天体重增加的预测概率是应用于我的特征(卡路里剩余量和步数)的加权线性组合的 sigmoid 函数加上常数截距项。这种关系在数学上表示为:
Superscript (i) denotes a a single data point/date
在哪里
- y(i): 日期 I 的体重再次状态(1 =体重增加，0 =体重减少)
- P(y(i)=1) :我在第一天体重增加的预测概率
- x(i): 在日期 i*时各个特征(卡路里剩余量和步数)的观察值
- θ :各个特征的回归系数/权重(从数据中学习)*
*请注意，我添加了一个额外的特征( x_intercept )，该特征将始终等于 1，以便可以学习截距项( θ_intercept )以及两个现有特征的系数。
因此，一旦从数据中学习到θ，就可以使用逻辑回归来分类我是否会在任何给定的一天增加体重——给定我的卡路里剩余量和步数——通过检查我在那一天增加体重的概率是否高于某个阈值(通常为 50%)。
如何学习 θ值
回归系数(θ)是通过最大化观察我的训练数据的概率的对数(也称为对数似然)来学习的。对数似然的公式是:
Log-likelihood of m training data points/dates
在哪里
- L :训练数据的对数似然(m 个数据点)
- y(i): 日期 I 的真实体重增加状态(1 =体重增加，0 =体重减少)
- P(y(i)=1) :我在第一天体重增加的预测概率(从之前的 sigmoid 函数中获得)
根据 sigmoid 函数，不同的θ值将产生不同的体重增加预测概率(P(y(i)=1)，从而产生不同的对数似然性。因此，我们的目标是找到最大化我的训练数据的对数似然性的一组θ，即最好地“解释”我的训练数据的θ。
使用批量梯度上升最大化对数似然
找到使我的训练数据的对数似然最大化的θ的一个简单算法是批量梯度上升，如下所述:
**第 0 步:**初始化θ的一些值(θ_ 截距，θ_ 剩余，θ_ 步长)
步骤 1: 对于每个训练数据点 I，使用该数据点的特征值(x_intercept*，x_surplus，x_step)和步骤 0 中初始化的θ的值，计算体重增加的概率。这是使用熟悉的 sigmoid 函数完成的:
*回忆 x_intercept = 1
步骤 2: 对于每个特征 j——截距/剩余/步长——使用下面的等式找到对数似然性相对于该特征θ的偏导数:
在哪里
- ∂L/∂θⱼ: 对数似然相对于特征 j 的θ的偏导数
- y(i): 日期 I 的真实体重再次状态(1 =体重增加，0 =体重减少)
- P(y(i)=1): 第一天体重增加的预测概率(来自步骤 1)
- xⱼ(i): 特性 j(截距/剩余/步长)在日期 I 的观测值，其中 x _ 截距(i) = 1，对于任何 I
步骤 3: 对于每个特征 j，通过对数似然性相对于θ的偏导数(来自步骤 2)乘以一个小常数(也称为学习率 α )来更新其θ:
这种学习速率控制算法收敛到对数似然最大值的速度，甚至控制算法是否收敛(有关更多详细信息，请参见后面的可视化算法收敛一节)。
步骤 4: 用这些更新的θ，重复步骤 1 到步骤 3，直到收敛。测试收敛的一种方法是查看对数似然是否已经收敛到稳定值，即它是否已经达到可能的最大值。
批量与随机梯度上升
步骤 2 中的求和符号——对所有数据点 I 求和 (y(i) - P(y(i)=1)) * xⱼ(i) 以计算偏导数——是该梯度上升算法属于批次种类的原因，因为每个偏导数是使用训练数据中的所有数据点计算的。没有这个求和符号，即偏导数仅使用一个数据点计算；假设数据点是随机选择的，梯度上升算法称为随机。
对于这个项目，我选择实现批量梯度上升，因为我的训练数据非常小(只有 46 个数据点)，所以一次使用整个训练数据集来计算偏导数没有问题。另一个原因是使用 numpy 的矢量化运算可以更容易地实现批量梯度上升(请参考下面的实现来了解如何实现)。
模型实现
数据预处理
Features & labels for 6 dates out of 46 in my training data
根据我之前的数据表，我使用剩余卡路里和步数列作为我的特征，体重增加状态列作为标签(见左)来训练我的逻辑回归分类器。
但是，在对该数据实施批量梯度上升算法之前，我首先需要:
- 在我的训练数据中添加一列 1 的来表示 x_intercept 特性的值。
- 通过(a)减去每个特征列的平均值和(b)除以其标准偏差，重新调整我的卡路里剩余量和步数特征。这样做有两个原因:
Sigmoid function
- 回想一下，sigmoid 函数(左)涉及 x 的线性组合的指数，因此 x 的非常小/大的值将使该指数内爆/爆炸。事实上，在我重新调整我的特征之前，我有可怕的趋同性，我不知道为什么。只有当我偶然看到 Jupyter 终端(而不是 Jupyter 笔记本输出)时，我才看到 numpy 默默产生的下溢/上溢警告页面！在我重新调整我的特征后，这些警告消失了，我的算法能够收敛了。
- 通过将我的特征缩小到相同的比例，我可以使用这些特征的回归系数(θ)来衡量它们在我体重增加/减少中的相对重要性。这将在我的模型的解释部分详细阐述。
在以上两个步骤之后，我的特征矩阵(X)被转换成维数为(46，3)的二维 numpy 数组，我的标签向量(y)被转换成维数为(46)的一维 numpy 数组。
Left: feature matrix X (46, 3). Columns (L to R): x_intercept, x_surplus, x_step. Right: label array y (46)
实施批量梯度上升
步骤 0:初始化θ的
theta = np.array([0.5, 0.5, 0.5])
我将所有θ初始化为 0.5。结果是一个维度为(3)的一维 numpy 数组theta
步骤 1:对于每个数据点，使用步骤 0 中的θ和 sigmoid 函数计算体重增加的概率
prob = 1 / (1 + np.exp(-X @ theta))
这就是 numpy 的矢量化运算派上用场的地方，因为 numpy 不是计算每个数据点的体重增加概率，而是一次计算所有数据点的体重增加概率:
X @ theta:通过将矩阵X与向量theta相乘，numpy 实质上是通过将X的每一行与theta列进行点积来计算每个数据点的 x 和θ的线性组合(见图中加粗的单元格)。1 / (1 + np.exp(-X @ theta)):在计算了所有数据点的 x 和θ的线性组合后，将 sigmoid 函数应用于它们中的每一个，以获得所有数据点的体重增加的最终概率。注意，这个 sigmoid 函数中的操作(/、+、-、np.exp)都表示 numpy 在后台运行的矢量化函数。这最终输出了概率向量prob，一个维数为(46)的一维 numpy 数组。
步骤 2:对于每个特征，使用步骤 1 中计算的概率，计算对数似然对相应θ的偏导数
gradient = (y - prob) @ X
这也是 numpy 的矢量化运算的亮点:
y - prob:对所有数据点的真实标签和预测概率进行直接的逐个元素的减法，得到表示差异的一维 numpy 数组(46 维)。(y - prob) @ X:在y - prob差向量(46)与特征矩阵X相乘之前，由 numpy 在幕后转置成一个维数为(1，46)的行向量，使其维数与X(46，3)的维数对齐。这种转置在 numpy 中也称为列向量的“向维度添加 1”。- 一旦尺寸对齐，
y - prob和X之间的向量矩阵乘法可以发生:点积在y - prob的行向量和X的每个特征列之间进行(见图中加粗的单元)。这实际上是对所有数据点 I 的**(y(I)-p(y(I)= 1))* xⱼ(i)**求和，以获得特征 j 的对数似然的偏导数(等式 2)。更令人印象深刻的是，在这种向量矩阵乘法下，可以同时计算所有三个特征的偏导数，从而产生一个称为gradient的偏导数向量。 - 从技术上讲，这个向量应该是(1，46)向量乘以(46，3)矩阵所得的维数(1，3)。但在幕后，numpy 在乘法后“去掉了(乘法前的)前置 1”，最后的
gradient向量是一个维数为(3)的一维数组。numpy 在乘法前后应用于其数组的这些幕后“扭曲”可以参考[numpy.matmul](https://docs.scipy.org/doc/numpy/reference/generated/numpy.matmul.html)的文档，它实现了矩阵乘法运算符@。
如果第一个参数是一维的，则通过在它的维数前加上 1，将它提升为矩阵。在矩阵乘法之后，前置的 1 被移除。
— numpy.matmul 文档
步骤 3:对于每个特征，用步骤 2 中的偏导数乘以学习率α来更新其θ
alpha = 0.01
theta = theta + alpha * gradient
这再简单不过了:我们可以通过将偏导数的gradient向量乘以某个预定义的学习速率alpha，并将乘积加到theta向量，来一次更新所有 3 个θ，而不是更新每个θ。
步骤 4:重复步骤 1 至 3，直到收敛
这可以通过将前面的 3 行代码嵌套到一个多次迭代的for循环中来轻松实现——每一步一行(见下面加粗的代码块)。下面展示的是运行 100 次迭代的批量梯度上升的完整算法(实现如此简单，难道不令人惊讶吗？！):
theta = np.array([0.5, 0.5, 0.5])
alpha = 0.01for _ in range(100):
**prob = 1 / (1 + np.exp(-X @ theta))
gradient = (y - prob) @ X
theta = theta + alpha * gradient**
检查收敛
检查算法收敛性的一种方法是查看对数似然性的差异在循环的过去几次迭代中是否保持在某个小容限水平以下，这表明对数似然性可能已经达到其最大值。使用y和prob(及其各自补码)的点积，对数似然可以简单地计算为:
log_likelihood = y @ np.log(prob) + (1 - y) @ np.log(1 - prob)
另一种(也许更有趣)方法是运行算法一定次数的迭代，并观察对数似然性是否达到了可能的最大值。下面是学习率α = 0.01 的批量梯度上升算法的 60 次迭代的可视化:
Left: convergence of θ’s. Middle: convergence of average log-likelihood. Right: convergence of decision boundary
- 从上面的动画中可以看出，平均对数似然性(对数似然性除以训练数据的数量，见中间的面板)在前 10 次迭代中快速上升，但之后开始趋于平稳。在第 60 次迭代时，迭代之间的平均对数似然的差异在 10^-5 的数量级，这表明足够好的收敛。这也对应于左图中回归系数(θ)的收敛。
- 另一种可视化这种收敛的方法是通过分类边界(见右图)。分类边界，也称为决策边界，代表体重增加的预测概率为 50%的线:其上的每个点的预测概率低于 50%，因此被分类为体重下降，其下的每个点的预测概率高于 50%，因此被分类为体重增加。在这 60 次迭代中，决策边界似乎稳定在一个合理的边界上，该边界将大多数体重减轻的日子(绿色)与体重增加的日子(红色)分开。
现在我们有了一个融合的逻辑回归模型，可以很好地对我的训练数据进行分类(至少通过目测决策边界)，让我们看看如何改进它。
模型改进
选择合适的学习速度
给定固定的迭代次数，学习率(α)的值可以确定算法在这些迭代之后是否会收敛，或者是否会收敛:
- 当α减少到 0.001(从最初的 0.01)时，θ的学习发生得更慢，并且在 60 次迭代之后，平均对数似然性仍然显示出增加的迹象。因此，在这个小的学习速率下，迭代次数应该增加，或者学习速率本身应该提高。
Dotted line represents the respective values at α = 0.01 for comparison
- 然而，如果学习率太高，θ可能会被“过度修正”，并在每次更新后在最佳点附近反弹。这一点从下面的α = 1 时的收敛动画可以看出，60 次迭代后似乎并没有收敛。因此，在这些情况下，应该降低学习率。
Dotted line represents the respective values at α = 0.01 for comparison
当α = 0.01 作为我们学习率的最佳点时，我们当然可以增加迭代的次数，以确保我们的模型很好地收敛到最大对数似然。事实上，在 1000 次迭代时，迭代之间的平均对数似然的差异实际上为零。
使用岭回归减少过度拟合
尽管我的逻辑回归模型已经收敛到我的训练数据的最大对数似然，但它可能会过度拟合训练数据，即它从数据中学习得有点太好了。因此，该模型可能会对我在 2018 年收集的过去数据进行很好的预测，但如果我用它来预测我在 2019 年的体重增加，可能会很糟糕。
减少逻辑回归过度拟合的一个解决方案是使用回归的 L2 正则化版本(也称为岭回归)，它通过一个由θ的平方组成的项减去原始对数似然:
m: number of training data points, n: number of features
因此，最大化上述函数相当于尽可能最大化训练数据的对数似然，同时保持θ较低(因为较高的θ将降低 L)。λ符号(λ)表示θ保持较低的程度(通常称为模型的正则化超参数)。当λ = 0 时，岭回归精确地返回到原始的、非正则化的逻辑回归。
在岭回归的实现方面，与原始逻辑回归的唯一区别在于批量梯度上升的偏导数(等式 2)的计算，其中从每个特征 j 的偏导数中减去 λ*θⱼ 正则化项:
需要注意的是，这种正则化并不经常用于截距θ，因此θ_intercept 的偏导数的计算与非正则化版本相同，即不减去λ*θ项。
该实现可以通过以下方式轻松集成到我们现有的 Python 代码中:
- 将
theta乘以正则化超参数lambda_reg得到reg_term——维数(3)的λ*θ正则化项。 - 将
reg_term的第一个元素设置为零，表示θ_intercept 没有被正则化 - 从
(y - prob) @ X中减去reg_term得到gradient。
以下是λ = 10 的岭回归代码，对原始算法的修改以粗体显示:
theta = np.array([0.5, 0.5, 0.5])
alpha = 0.01
lambda_reg = 10for _ in range(100):
prob = 1 / (1 + np.exp(-X @ theta))
**reg_term = lambda_reg * theta
reg_term[0] = 0**
gradient = (y - prob) @ X **- reg_term**
theta = theta + alpha * gradient
当在 60 次迭代中监测该岭回归的收敛时，我们可以看到:
Dotted line represents the respective values at λ = 0 (non-regularized logistic regression) for comparison
- 在λ = 10 时，我的两个主要特征(卡路里剩余量和步数，见左图)的θ收敛到更接近于零的值，即与非正则化版本(λ = 0)相比幅度更低。然而，截距θ的收敛基本不受影响。
- 岭回归的平均对数似然比非正则化版本收敛到更低的值(见中间面板)，这表明岭回归对我的训练数据提供了不太完美的拟合，但这也可能意味着它对我的训练数据的过度拟合更少。
- 岭回归的决策边界稍微偏离非正则化边界(见右图)。然而，它似乎仍然能够很好地区分我的训练数据点(红色和绿色)。
随着λ的增加，学习回归系数(θ)进一步向零压缩，截距除外(见下图)。此外，从决策边界可以看出，岭回归在对我的训练数据进行分类时变得越来越不有效:λ = 1 时的决策边界非常接近非正则化边界，而在λ = 100 时，边界实际上不可用(见右图)。
How θ’s and decision boundary changes at different λ values
然而，岭回归的目的不是改善对训练数据的拟合(因为如果是这样，它的表现将总是比非正则化版本差，如上所述)。更确切地说，它是为了改进对新数据的预测，即没有训练过岭回归的数据。
为了比较岭回归和它的非正则化对应物，我使用了双重交叉验证，如下所示:
1。 将46 个数据点随机分成 2 等份:A & B(各 23 个数据点)
2。 在部分 A——训练集——上训练岭回归，并在训练集上记录回忆(正确预测体重增加的天数/部分 A 中真实体重增加的天数)
3。使用在 A 部分训练的θ来预测我是否会在 B 部分(验证集)增加体重，并在验证集上记录回忆(正确预测体重增加的天数/B 部分中真实体重增加的天数)
4。 交换零件，重复步骤 2 和 3，即零件 B 现在是训练集，零件 A 是验证集
5。 对两次试验中的训练集召回进行平均，验证集召回也是如此
为什么要用回忆？
有两个基本的度量标准来衡量分类器工作的好坏:精确度和召回率。在这方面:
- 就我个人而言，如果分类器预测我会增加体重，而事实证明我不会，我不会在意；事实上，这甚至会是一个受欢迎的惊喜！所以在预测的增重天数中，我并不是太在意假阳性。换句话说，精度对我来说没那么重要。
- 另一方面，我更可能纠结于分类器是否会标记出我真正增重的所有日子，以免它们被误报为假阴性(分类器预测我会减肥，而事实恰恰相反)。换句话说，我会尽可能提高召回率。
以下是不同正则化水平下的平均训练和验证召回，λ范围从 0(非正则化)到 10:
Average train & validation set recall (across 2 folds) at different λ values
- 从左图来看，随着λ从 0 开始增加，平均训练集召回率保持在恒定水平，并且仅在λ接近 10 时开始下降。这与早期的观察结果一致，即岭回归对训练数据的性能随着λ的增加而变差。
- 另一方面，对于λ在 0.1 和 1 之间，平均验证集召回具有明显的隆起，这表明在这些λ处的岭回归在新的验证数据上比其非正则化对应物表现得更好，即使两者在训练数据上表现得一样好。
- 通过进一步的检查，发现在一个验证组中，真实的体重增加日(红色带黑边，见下图)被非正则化回归错误地归类为体重减少:它停留在(实)决策边界之上。另一方面，当λ = 0.5(介于 0.1 和 1 之间)时，决策边界会略微扭曲。因此，该点保持在该(虚线)决策边界之下，并且被正确分类。这是岭回归在λ = 0.5 时比其非正则化对应物在验证集上具有更好召回率的唯一原因。
Left: Decision boundary on one validation set (with training set faded in background). Right: dotted region on left panel enlarged
尽管这些结果表明我应该选择具有最高平均验证集召回率(λ = 0.5)的λ值，但验证集(23)中的少量数据点(其中甚至有更少数量的重量增加点(红色))表明这种性能的提高可能只是由于运气。这也解释了为什么对于一些λ，在相同的λ下，验证集的性能高于训练集的性能，即使通常发生相反的情况；我可以有一个“幸运的”验证集。
也就是说，选择λ = 0.5 没有害处，因为当在我的 46 个点的整个数据上训练时，它的决策边界实际上与λ = 0 的决策边界没有什么区别，正如前面不同λ(从 0 到 100)的决策边界图所示。因此，对于我的最终模型，我选择保持λ = 0.5。
选择决策边界的阈值
在调整我的模型时，选择λ并不是最有影响的决定。相反，选择我的决策边界的阈值是:
- 对于我到目前为止建立的所有回归模型，分类阈值，也就是决策边界，被设置为 50%(或 0.5)。这是一个合理的阈值，因为预测体重增加概率高于 50%的日子自然应被归类为体重增加的日子，反之亦然。
- 然而，当这个体重增加阈值降低时，越来越多的数据点将被归类为体重增加(正确与否)。回忆，也被称为真阳性率，当然会增加，但是假阳性率——在真正体重减轻的日子里被错误地预测为体重增加的日子——会同步增加(见左图，以及下面中间图中的 ROC 曲线)。然而，如前所述，由于我不太担心假阳性，我可以容忍许多体重减轻日被错误地归类为体重增加，即高假阳性率，如果这意味着我的真实体重增加日被更好地检测，即高真阳性率/回忆。
Left: True/False Positive Rate at different classification threshold. Middle: ROC curve. Right: decision boundary at different thresholds
如果是这样的话，门槛是不是应该定得越低越好，甚至到 0%？不，因为:
当阈值降低时，决策边界上移。
- 例如，为了使回忆增加到最近的较高水平，阈值必须从 50%减少到 44%(在左侧面板中从黑色减少到棕色)。结果，决策边界上移，以捕捉更多的体重增加点(右侧面板中带有棕色边框的红点)。这相当于向左移动了 107 千卡。换句话说，如果我以前在 50%的界限，我必须在相同的步数下少吃 107 千卡，才能保持在 44%界限的良好减肥侧。
- 在阈值接近 0%的极端情况下(左侧面板中的橙色点)，边界向上移动，以便捕捉所有体重增加点，包括带有橙色边框的最高红色点(右侧面板)。这个决策界限决定了在我懒惰的一天，我平均走了将近 2500 步，我应该比我的预算少吃 1740 千卡来保持健康。鉴于我的数据中的平均预算约为 1715 千卡，这转化为在那些懒惰的日子里的上限 -25 卡路里(是的，你没看错)。当然，在这个物理极限下，我的回忆会是惊人的，但是我会死的！
- 那么我到底该不该降低自己的门槛呢？鉴于从我的数据中学到的 50%决策界限已经相当保守了——在懒惰的日子里，这表明我应该比通常的 1700+千卡预算少吃大约 140 千卡——我决定坚持默认的 50%阈值。将这个阈值降低到下一个可能的水平 44%会让我在回忆中多得几分，但额外的 107 千卡限制不值得我失去理智。
模型解释
概括地说，我们通过梯度上升算法训练了一个逻辑回归分类器，使用我的每日卡路里盈余和步数作为特征，我的体重增加状态作为标签。我们模型的参数是:
- 学习率: α = 0.01
- 正则化参数: λ = 0.5
- 阈值: 50%
在整个数据集上训练分类器后，学习的回归系数(θ)得出:
优势比分析
- 逻辑回归系数的一个常见解释是通过 优势比 :特征的一个单位变化，结果变化的几率是多少倍。对于一个特征 j，那个特征的比值比正好是它的θ的指数。
- 回想一下，逻辑回归是对卡路里盈余的标准化值进行的(通过减去盈余平均值并除以盈余标准偏差得到)。因此，θ_surplus = 0.9 意味着在任何步数下，剩余卡路里的标准差减少(约 420 千卡)对应于 e^0.9，或者我体重增加的几率减少 2.5 倍。
- 另一方面，在θ_step = -1.2 的情况下，我步数的标准差增加——大约 5980 步——对应于 e^1.2，或者我体重增加的几率减少 3.3 倍(在任何卡路里剩余量的情况下)。在我最近的 10 次锻炼中，我的平均步频是 1366 步/公里。因此，这 5890 步相当于大约 4.4 公里。
- 换句话说:
步数的标准差增加(5980 步)比卡路里消耗量的标准差减少(420 千卡)更能有效减少我体重增加的几率。
当然，我愿意跑 4.4 公里而不是不吃那碗 420 千卡的河粉完全是另一回事！
尽管比值比为我应该如何制定减肥策略提供了一些见解，但从我的逻辑回归分类器的决策边界中可以得出一个更可行的计划。
决策边界分析
回想一下我们之前可靠的 sigmoid 函数:
Superscript (i) removed for simplicity
很容易看出，体重减轻的概率为 50%(其中分类阈值为)，θ和 x 的线性组合必须为零*(我还用 1 代替了 x_intercept):
50% decision boundary formula for normalized features
*对于其他概率阈值，θ和 x 的线性组合可以通过取阈值的 logit 得到:ln(p/(1-p))
然而，在我们的逻辑回归中使用的 x 是我们原始特征的标准化值。因此，我们可以将上面的等式改写为:
***** denotes the original feature value. μ: feature mean, and σ: feature standard deviation
重新排列，我们有:
50% decision boundary formula for original features
将θ和特征均值(μ)和标准偏差(σ)代入上述方程，得到原始卡路里剩余量和步数特征之间的线性方程:
从图形上看，这个方程代表了原始剩余步骤图中的决策边界(见下面的左图)。根据这个决策界限，我应该记住 3 个重要数字(见右图):
Left: surplus-step plot with decision boundary. Right: dotted region from left panel enlarged
1) -140 千卡
这是我在懒惰的日子里(平均 2480 步)应该吃的低于 Loseit 应用程序卡路里预算 的量。平均预算刚刚超过 1700 千卡，这意味着在那些日子里，我应该平均摄入低于 1560 千卡的热量。这个数字听起来确实很严格。****
130 千卡
然而，一个可取之处是，根据决策边界，对于我超出正常活动范围的任何 1 公里，我都可以将这个限制增加 130 千卡。例如，如果我计划在某一天跑 5000 米，我可以负担得起比应用程序当天的卡路里预算多吃 140 + 130 * 5 = 510 千卡。希望这能鼓励我坚持我的跑步计划。
3) 1070 步
另一方面，当我试图吃得超过卡路里限制时(如前两条规则所规定的)，我吃的任何超过限制的 100 千卡必须通过至少 1070 步来获得。这可以通过以下方式实现:
Hit me with that number baby!
结论
我希望用以上的指导方针，让我 2019 年的减肥之旅比 2018 年更成功。当然，即使有数字支持，我也不会总是成功，但我从这个项目中学到的最重要的一课是:
我应该善待自己。
例如，我减肥，即使我在活跃的日子里吃得比预算多，所以我不应该为这样做感到内疚。希望有了这种结合饮食和锻炼的新模式，在我持续的减肥旅程中，我可以感觉到更少的罪恶感和羞耻感(就像我以前一样)。
其他课程
我着手这个项目的另一个重要原因是实现一个机器项目，而不使用预先存在的库，如 scikit-learn。以下是我从这个过程中学到的一些经验:
- 一年多前，当我第一次学习 Python(或者任何严肃的编程)时，我不确定为什么会有人使用类。嗯，在这个项目中(a)使用全局变量
theta, (b)得到一些奇怪的结果，比如不收敛，以及©意识到所述的theta属于我许多个月前运行的其他模型，我现在意识到为什么类中的封装是如此重要:一个模型对象可以有它自己的属性(theta、alpha、lambda_reg)和方法(fit、predict)，这些属性和方法不会与其他模型对象的属性和方法冲突，我可以在下次需要它的时候愉快地把它完整地捡起来。
我真的对如何使用类来建模数据科学问题很感兴趣，并且我认为通过更多相关的例子，我可能能够更好地理解面向对象编程的强大功能(现在我仍然不太确定我应该何时或者如何使用它们)。
—初学 Python 时的我
- 这个项目还允许我实现和理解一些我从未有机会使用的简洁编程概念的实际好处，例如使用生成器即
yield语句一次返回一个训练和验证文件夹。这允许我仅在这些折叠上评估我的模型，然后返回并生成更多的折叠，而不是一次将所有折叠存储在内存中(尽管对于我的数据大小来说，这几乎没有区别)。甚至像编写清晰的伪代码这样的事情也被证明是非常重要的，尤其是在实现像逻辑回归这样的数学算法时；令人尴尬的是，我不得不承认我花了很长时间进行调试，因为我没有事先编写明确的伪代码，并且错误地交换了alpha和theta。 - 最后，这个项目帮助我更详细地研究了我经常使用的一些机器学习库的实现，例如偶尔警告我在训练 scikit-learn 模型时必须指定
max_iter或tol:前者指定让模型收敛的迭代次数，后者指定迭代停止的容差水平——这正是我在为自己的模型检查收敛性时面临的两个选择。 - 另一个例子是，我现在了解到一些 scikit-learn 模型的行为不像我期望的那样:
SGDClassifier(loss=’log’, penalty=’l2', learning_rate=’constant’)似乎没有缩小截距的θ，并给了我的模型类似的θ，而LogisticRegression(penalty=’l2')在默认设置下缩小截距的θ，除非有人篡改了intercept_scaling参数。因此，我在将来使用第三方库分析我的数据时会更加小心，并在必要时验证结果。
更多资源
- 我从 CS229 的讲义和视频 讲座(我极力推荐的吴恩达教授讲授的机器学习课程)中获取了逻辑回归及其梯度上升法的大部分推导。注释提供了随机梯度上升的公式，但是批量版本可以很容易地修改(如前面的数学回顾部分所示)。
- 一个更容易理解的解释可以在他在 Coursera 上的第三周机器课程中找到。本课程也涵盖岭回归在逻辑回归上的应用，这是 CS229 课程所没有的。一个小细节:Coursera 课程将逻辑回归称为使用随机梯度下降最小化对数损失。然而，这与我实现的使用随机梯度上升最大化对数似然性是一样的，对数损失只不过是对数似然性的负值(有一些小的修改)。
- 对逻辑回归和岭回归的另一个很好的解释来自华盛顿大学 Coursera 关于分类方法的课程。幸运的是，本课程使用对数似然最大化来解释逻辑回归，所以它应该与 CS229 和我的符号一致。
- 最后但并非最不重要的是，我在网上找到的两个减肥分析项目真正启发了我自己的项目:一个来自威尔·科尔森，他只使用他过去的体重来预测未来的减肥，另一个来自阿里尔·费根，他使用几十个因素来预测他每天增加或减少多少体重，并得出非常有趣的结果:显然睡眠是他减肥的最重要因素！我的方法介于两者之间，仅使用两个特征来预测我是否会减肥，这使得我可以在我的减肥旅程中带走简单的可视化和简单可行的见解
我希望我的项目可以启发其他人使用机器学习和数据科学来帮助他们更多地了解自己，完成减肥等个人目标。如果您有任何问题或反馈，请不要犹豫通过媒体联系我！
分析网易音乐-第一部分:播放列表
中国人听什么样的音乐？
Netease Music Logo
网易音乐(https://music.163.com/)，一个中国版的 Spotify，是一个音乐应用程序，允许用户播放、下载和参与来自世界各地的音乐。网易的月活跃用户数超过 7000 万，与酷狗音乐、酷我音乐和 QQ 音乐并列中国音乐应用前四。此外，它是四款应用中唯一一款非腾讯的应用。
网易音乐的受众呈现出一些有趣的特征。根据该公司的报告以及 QuestMobile 的数据，该应用的用户更年轻，男性比例更高，并经常与社区互动。用户不仅比其竞争对手更频繁地打开和使用该应用，而且他们也倾向于留下大量评论。事实上，网易音乐中的评论已经成为其一大看点。人们会发布有趣的笑话、心碎的故事，甚至向他们的秘密情人表达自己的感受，并要求其他听众竖起大拇指，这样他/她就可以看到了。在网易音乐的一次推广活动中，它在杭州的地铁站和火车上发表了点赞数最高的评论。这一创意活动是该应用程序的巨大成功，并进一步确立了其作为最具“社交”属性的音乐应用程序的地位。
The comment promotional campaign
与西方同行类似，大多数听众将使用播放列表作为存储他们喜爱的歌曲以及发现新事物的方式。在网易音乐中，用户可以“播放”该列表，“收藏”该列表以供将来收听，“评论”该播放列表，并在应用程序内或微信/微博中“分享”该播放列表。在本文中，我将单独查看播放列表的一些指标。在接下来的几篇文章中，我还将探索关于单首歌曲的指标以及这些歌曲的评论。
如何获取数据
所有的数据都可以在网站上公开获得，并使用我个人构建的爬虫收集。爬虫基本遵循从主页发现播放列表，从播放列表发现歌曲，然后从播放列表和歌曲下载评论的逻辑。爬虫仍在积极开发中。如果你对学习如何爬行感兴趣，请在下面留下评论，我会写一篇关于爬行的文章。
数据大小
下面进行的分析使用了 6000 个播放列表的样本。这些列表要么是在主页上推荐的，要么与那些热门播放列表有些关联，所以这些列表往往有更多的参与。
探索数字变量
首先，我们将查看数据集的数值。如前所述，参与度指标是“收藏”、“评论”、“分享”和“播放”。数据集还有一个“song_num”变量，它描述了播放列表中歌曲的数量。由于默认的 Pandas describe()使用 float 类型并在表中留下许多 0，我们可以使用 applymap()函数更改输出的格式，使它更漂亮。
df[['comment_num', 'fav_num', 'share_num', 'play_num','song_num']]
.describe()
.applymap(lambda x: format(x, '.0f'))
describe() result
从这里我们可以看到，评论和分享的中位数是 10。“fav”的号码是 663，播放号码是 36k。粗略地说，一个播放列表需要播放 4000 次才能获得“分享”或“评论”。这表明评论和分享功能并不常用。请注意，这里的“评论”适用于整个播放列表，它不同于评论一首更受欢迎的歌曲。
为了直观地显示这些指标，我们可以使用一个方框来绘制这些列。然而，由于所有四列的数据都有很大的偏差，直接绘制它们不会产生好的结果。下图并没有真正告诉我们什么，因为 x 轴非常大，没有显示足够多的“长尾”播放列表的细节。
sns.boxplot(data = df[['comment_num', 'fav_num', 'share_num', 'play_num']], orient='v')
Boxplot without normalization
为了使这种偏斜数据正常化，我们可以使用对数函数来做一个快速而肮脏的把戏。如果你想了解更多细节，请阅读这篇文章。
https://becoming human . ai/how-to-deal-with-skewed-dataset-in-machine-learning-AFD 2928011 cc
sns.boxplot(
data = df[['comment_num', 'fav_num', 'share_num', 'play_num']]
.applymap(lambda x: math.log(x+1)), orient='v')
Boxplot with normalization
归一化后，我们可以看到所有四个指标都遵循相对正态分布。
变量之间的相关性
播放次数多的播放列表也有很多评论吗？这似乎是一个明显的正相关，但看看这个数字如何滚动出来仍然很有趣。
corr = df[['comment_num', 'fav_num', 'share_num', 'play_num']].corr()
sns.heatmap(corr,
xticklabels=corr.columns,
yticklabels=corr.columns,
cmap = 'Blues',
annot = True)
我们可以使用 pandas 的 corr()函数来获得相关图，并使用 seaborn 来绘制相关图。注意，您可以使用 cmap 参数来改变图形的颜色。对于关联热图，基于单一颜色的图表会更有效，因为它更清楚地显示了数字比较。
有意思的是，播放号和收藏号的相关性最强。这是有意义的，因为在他们可以保存播放列表之前，需要有大量的人在听它，收藏它意味着他们将来会听它。评论与其他指标的相关性较弱。由于评论数和分享数与其他两个相比要小得多，很可能这些指标中有较高比例的“噪声”,因此导致相关性较弱。
另一方面，股票数量的相关性更弱。尽管仍有很强的正相关性，但这个数字并不那么确定。也许一个原因是，共享播放列表的意义不同于播放它或喜欢它。通过演奏或收藏，听众可以独立地欣赏音乐。通过对播放列表进行评论，听众可以与作者和喜欢该播放列表的其他人进行互动。然而，通过共享播放列表，这意味着听众在社交上参与音乐，并致力于向他们的朋友推荐该列表。也许一个人会在白天分享一个经典的摇滚播放列表，然后在晚上循环播放“我知道你是个麻烦”?
探索其他变量
我们可以看的另一个有趣的列是“tag_list”。当用户创建播放列表时，用户最多可以选择 3 个标签。这些标签包括语言、流派、合适的场合、情感和主题。对于那些不懂中文的人，我创建了一个翻译文件来将这些标签翻译成英文。
Tag selection page at Netease Music App
在将标签列表展开成一个带标签的标签列表后，我们可以绘制每个标签的直方图。虽然样本量很小，但数据也许能显示出网易音乐听众音乐偏好的一些有趣趋势。
语言
Language count plot
在语言方面，中国音乐仍然是最主要的类别，但西方音乐(主要是英语)紧随其后。日语也占据了播放列表的相当一部分。
流派
Genre count plot
就流派而言，除了很难精确定义的“世界音乐”之外，流行音乐是最受欢迎的一类，因为人们不知道还能选择什么。“新时代”和“轻音乐”都很放松，占据了第三和第四位。这可能会令人惊讶，因为这些音乐在西方世界并不流行。尽管 EDM 和 Rap 最近在年轻一代中很受欢迎，但它们在播放列表中所占的比例仍然相对较小。
场合
Occasion count plot
对于场合，夜晚和午休是两个占主导地位的标签。这与我们上面确定的“轻音乐”趋势相一致，因为人们在这两个时间段倾向于听到更多的“舒缓”音乐。同样值得注意的是，对于“流行”列表，列表创建者可能不会为其选择“场合”标签。
情感
Emotion count plot
Theme count plot
排在前两位的情绪是“怀旧”和“开心”。与“悲伤”和“孤独”列表的低数量相比，人们似乎更倾向于制作积极情绪的播放列表。就主题而言，由于数量很少，很难确定一个趋势。值得注意的是，动画和游戏具有很强的代表性，表明“动漫”文化在中国仍然蓬勃发展。
总体而言，与西方世界相比，网易音乐听众的口味更“清淡”。观众作为一个整体已经接受了一套不同的音乐。
有趣的事实
- 播放次数最多的播放列表播放次数达 1.48 亿，是一个收录了最佳 intros 的英文歌曲集。榜单第一首是《五百里》。
- 评论最多(10k)的播放列表是 90 年代中国动画片的配乐集。
好了，这一次到此为止。下一篇文章，我会在播放列表数据上探索一些机器学习的方法，尝试基于文本数据和歌曲数据来预测类别。
分析曼哈顿一百万选民的记录
对曼哈顿选民构成的独特观察
如果你能立即想象整个城市的政治归属，下至每一个公寓和登记投票的人，会怎么样？有点令人惊讶的是，纽约市在 2019 年初将这一点变成了现实，当时纽约市选举委员会决定在网上发布460 万选民记录，由 《纽约时报》 报道。这些记录包括全名，家庭住址，政治派别，以及你在过去 2 年中是否登记过。根据此条的原因是:
委员会官员表示，他们的印刷供应商无法及时制作足够的选民登记册，供候选人在 2 月份开始收集请愿书签名，这就是他们在网上发布信息的原因。
所以这成了其间的‘权宜之计’。尽管正如纽约市 BOE 发言人瓦莱丽·瓦兹奎·迪亚兹在同一篇文章中所说:
不管你有没有意识到，纽约的选民数据是公开的。
虽然一个城市心血来潮公开 460 万份选民个人记录的想法令人深感不安，我打算在另一篇文章中讨论这一点，但可视化这些数据的机会引起了我的兴趣。对于大多数公共数据集，你可能能够访问的最细的数据集可能是人口普查区级别，但是纽约市刚刚展示了世界上最独特的数据集之一！现在第一次有了一个数据集来可视化大量超粒度选民数据，包括按地区划分的政治派别集中度，在某些地区是否有协调一致的努力来登记新选民，以及占主导地位的政党在多大程度上是人们可以理解的。
最终我决定尝试并发表这篇文章，原因有二:
1。公共官员需要更好地考虑在线发布这样一个大规模数据集的潜在后果
在这种情况下，一个半技术人员花了几个周末的空闲时间来下载、转换和搜索数百万条记录。我不仅做到了这一点，而且在谷歌上获得积分后，我能够以大约 30 美元的价格对超过 100 万条记录进行地理编码。在这个数字时代，隐私需要放在首位，而不是事后才想到。
②。总的来说，这些数据提供了新的、独特的见解，对很多人都很有价值**
我确实相信这个数据集(在建筑层面可视化/分析)可以提供独特的见解，同时维护我自己和其他纽约居民的隐私。
总的来说，我发现确实存在政党主导的集中区域，而且政党主导的集群遍布整个城市。我还发现，新选民登记高度集中在字母城和上东区等地区，一些集群的登记人数比其他集群高 5 倍。本文的其余部分将涵盖将大约 600 个 pdf 转换成可用数据集的方法、观察到的独特比较的三种可视化，以及对进一步研究提出建议的结论。
方法学
现在，仅仅因为数据可用并不意味着它是可访问的，数据发布的方式证明了这一点。让我们来看看下面经过编辑的原始 PDF:
Redacted sample of the PDF pages
从数据科学的角度来看，这绝对是一场噩梦。许多原因中的几个:
- 它在 PDF 中显示为每页三个“表格”,因此 excel/其他程序不会像您预期的那样使用它。它或多或少会爆炸，并在随机单元格中生成随机条目，因此不适合简单的复制/粘贴
- 数据重叠(街道名称与名称重叠),因此没有简单的方法来分隔列
- 像街道号这样的数据点最初只陈述一次
- 列/新页面可以从没有数据开始(例如，第二列没有街道号码),因此需要有一种方法来记住“最近”的元素
- 它只有地址，为了在地图上可视化，您需要纬度/经度数据
由于我是 Python 爱好者，我立即开始尝试使用各种现有的库，如 PyPDF2 来尝试读取数据，但是很明显 PDF 的格式不允许这样做。在陷入死胡同的最初几个小时后，我陷入了困境，开始重新思考我的策略。如果我不能使用页面上现有的数据结构，我需要找到一种方法来完全隔离和重建页面上的每个单词。
输入光学字符识别(“OCR”)。使用 Pytesseract 库，我改变了策略，告诉我的程序每一页都是图像而不是 PDF 页面。事实上，只需要一行代码:
text = pytesseract.image_to_data(第页)
为您提供页面上几乎每个元素的坐标:
Output of Pytesseract using image_to_data
一旦我确认了这一点，接下来的问题就是如何使用这些位置来构造一个自定义的分隔符来遍历每一页。不幸的是，由于 Pytesseract 选择的分组方式，前几列中的大部分都不起作用，但据我所知，“左”、“上”、“宽”和“高”的信息是准确的。所以对于我的程序，我基本上利用了以下逻辑:
- 将页面作为图像阅读
- 隔离左侧的数据表
- 使用左/上位置来确定行/列分组
- 制作一个数据框架，其中包含街道号、姓名、邮政编码等单独的列。
- 对中间和右侧的数据表执行上述步骤
- 移到下一页
公平地说，我不得不处理许多细微差别和边缘情况，OCR 的结果并不完美(我将在结论中详细讨论)。然而，经过一两天的修补，我能够在大约 30 分钟内成功地迭代通过我的第一个 250 页的 PDF，并将其转换成一个完全可利用的数据表！
Sample view of our newly converted dataset
现在是有趣的部分，我如何把它翻译成实际的地图？我有地址，但是我不能实际映射它们，除非我能以某种方式找到每个地址的纬度/经度坐标。当我查看我成功转换的唯一曼哈顿地址时，它达到了近 100，000 个唯一地址，因此手动完成这一部分是不可能的。幸运的是，谷歌不仅有一个反向地理编码 API，而且我发现了一个名为 Geopy 的奇妙的 python 库，它让我可以轻松地创建一个函数，用不到 20 行代码对地址列表进行地理编码。
About 20 lines of code can geocode nearly any address
由于我从未使用过谷歌的云平台，我从 300 美元的免费信用开始，这有助于抵消大约 330 美元的总成本。在我的免费积分和奇思妙想的鼓舞下，我让我的程序在接下来的 12 个小时内执行，为几乎所有相应的地址返回一组漂亮的经纬度坐标。在将它与聚合数据集连接起来并进行一些数据清理之后(见结论)，是时候决定如何可视化它了。
可视化大规模数据集不是一件容易的事情，除了找到一个可以处理几十万行的可视化库之外，您还需要考虑可视化数据集的技术，以便最终结果对底层数据来说是有意义和真实的。我最初开始在机器学习社区使用一个奇妙的新项目，名为 Streamlit ，但最终决定我不想创建一个仪表板，因为太多的杂物会分散注意力(绘制数据)。
First Attempt with Streamlit, Great for Creating Dashboards for Active Analysis!
最终，我选择了优步可视化团队制作的名为 Kepler.gl 的图形库。开普勒使用 Mapbox 作为底层地图提供商，不仅在处理大型数据集方面做了惊人的工作，而且即使对于不编码的人也非常容易访问。在尝试了一些不同的技术、滤镜和图层设计后，我决定选择下面三个最有影响力的。
选民的绝对优势
我想看的第一件事是’孤立地看每个聚会浓度是什么样的’？我想这样做有两个原因:
- 为了看看我提取的数据是否大致接近报道的党派优势(民主党已经占了大约 68%的投票人口)
- 为了清楚地识别集中区域，没有一方由于偏斜而控制地图
为了做到这一点，我在地图上创建了每个政党的双重视图，根据给定区域中的点的集中程度(该区域中注册选民的密度)创建了半径宽度，然后根据 avg 创建了连续的量化色标。双方所占面积中的建筑百分比。
Democratic (left) concentrations compared to Republican (right)
正如我们在地图对比中看到的那样，曼哈顿确实是压倒性的民主，许多倍数更集中在像哈莱姆区(曼哈顿北部)这样的地区。毫不奇怪，他们往往也是大多数建筑中的大多数选民，正如你在地图上看到的深蓝色。从共和党方面来看，很明显，从原始数据来看，他们占少数，而在他们拥有大量集群的地区，民主党人通常仍占主导地位(由浅红色和深红色显示)。与其他地区相比，上东部地区，如鲁诺山、海龟湾和约克维尔确实显得高度集中。林肯广场和金融区周围的人群似乎也较少，尽管要稀疏得多。
选民的混合优势
现在我们已经检查了原始的政党优势，让我们试着合并到一张地图上。这里的想法是试图显示政党集中的主导领域，同时也保持人口的相对规模。为了做到这一点，我创建了一个网格样式的层，使用了 0.2 公里(约 0.06 英里)半径内的聚会总数。然后我会为每个网格创建一种颜色，基于在该区域找到的最频繁的聚会。我的假设是，即使民主党控制了曼哈顿，也可能有共和党控制的地区。下面是结果:
Party dominance in clustered areas
正如我们所看到的，这并不令人惊讶，因为民主党主导了地图，即使是在上东区这样的共和党人最集中的地区。似乎确实有一些主题区域，如中央公园的东南区域，但最终在同等规模的民主集中度下，几乎没有优势。
新选民登记
我特别好奇的一件事是，曼哈顿的选民登记是否集中在某些地区。pdf 包括 2017 年和 2018 年登记的选民数据，所以我决定看看可以收集到哪些见解。总的来说，我能够捕获大约 55，000 名新注册的选民，所以我决定将新选民聚集在 0.2 公里(~0.12 英里)的集群中。曼哈顿本身只有大约 2.3 英里宽，所以选择一个对你所代表的区域有意义的比例尺是很重要的。起初，我试图像上面的双地图一样将政党分开，但令人惊讶的是，在新民主党人和共和党人登记的地区，我看到的变化很小。相反，我决定将两党的新登记者结合在一起，并从总体上看 2017 年与 2018 年的新选民。
2017 New Voters (left) compared to 2018 New Voters (right)
正如我们所看到的，有两个新选民登记非常集中的地区——上东区的鲁诺山/约克维尔地区，东村的字母城周围地区。很难从中收集到任何明确的见解，但这种集中很有趣，因为地图上红色区域的新选民集中程度几乎是绿色区域的 5 倍。有几种方法可以潜在地解释这一点——也许这些地区相对于其他地区来说，由更多的新移民组成；或者，这可能反映了社区集中努力鼓励选民登记。
结论
这是一个非常有趣的项目，但是我还需要强调一些关于分析的免责声明。首先也是最重要的，像我这样使用 OCR 在准确性方面有很大的折衷。我的最终数据集被精简为大约 750，000 条聚合记录，其中大部分是由于 OCR 读取错误而不得不丢弃的数据(例如/“100 Orange Street”可能会读作“T0O 0range $treet”)。
因为像街道号码和名称这样的数据在某些情况下只出现一次，程序以这种方式误读它们最终会导致纬度/经度坐标找不到或不正确。我尽我最大的判断，试图扔掉明显不正确的大量数据，但我知道这个数据集并不完美，可能无法 100%代表现实。虽然数据可能是准确的，但我没有完美的验证方法，因为这是一个人的操作。
我很乐意看到城市思考如何创建像这样更开放的数据集。虽然我对纽约市发布数据集的方式感到失望，但我认为汇总层面的基础数据对许多相关方来说都有巨大的价值。这不仅有助于研究人员更好地了解选民趋势，还可能极大地促进与选区划分和选区划分等问题相关的研究。
进一步的研究可以纳入其他数据集，如美国人口普查数据，以了解诸如地区选民参与度或刺激选民登记的社会经济混合等情况。我也没有想象像布鲁克林或共和党或民主党以外的其他政党这样的地区，所以这些地区也仍然是未开发的领域。我相信，深入研究这些领域会发现更多关于纽约选民构成的独特发现。我希望这能有所帮助，并给每个人提供思考的食粮——感谢您花时间通读，并在下面自由发表您的想法或问题！
分析在线活动和睡眠模式
让数据科学变得有趣
分析你的脸书朋友的网上活动和睡眠模式
社交网络上有大量公开的信息，有时我们甚至会忘记它们的存在。像我们的脸书朋友在的在线活动这样的小信息能让我们推断出他们什么时候睡觉，或者他们什么时候在脸书最活跃。
在本帖中，我们将学习分析脸书好友的在线活动，并找出使用统计和睡眠模式。该员额由两部分组成。
- 展示在线活动和睡眠模式
- 使用 Python 分析在线活动和睡眠模式(带代码)
在第一部分，我将向你展示一些我用我脸书朋友的数据制作的有趣的图表，在第二部分，我将向你解释如何用 Python 创建一些展示的图表。
本帖完整代码的 Jupyter 笔记本可以在 这里 找到。
1.展示在线活动和睡眠模式
从一个人的脸书朋友的在线状态数据中可以得出一系列有趣的见解。在这一节中，我将向你展示一些有趣的见解，这些见解来自我通过爬行脸书收集的两周数据。
首先，我非常有兴趣找出我的朋友们在一天中的不同时间以及一周中的不同日子使用脸书/Messenger 的模式。
给定一个时间戳列表和相应的当时在线的朋友列表，我们不仅可以找出一天中最多和最少数量的朋友在线的高峰时间，还可以对 睡眠模式 有一个相当好的估计。
这些统计数据也提醒我们，我们在网上留下了多少公共信息和痕迹，以及它们是如何被以不同方式利用的。
虽然数据收集非常简单，但是与收集的数据样本及其预处理相关的一些要点值得一提，如下所示:
- 在 两周 的过程中，每周 7 天，每天 24 小时，每 4-5 分钟收集一次样本。仅仅两个星期的短时间跨度可能导致有偏差的数据集，但是正如我们将在后面看到的，它给了我们一个很好的模式估计。
- 我的大多数朋友都是穆斯林，数据也是在神圣的斋月收集的，在不同的时区(包括欧洲、美国、亚洲和非洲)都有朋友。这可能会导致睡眠模式与平时不同。
- 当比较朋友的相对在线活动时，使用了每个朋友的当地时区例如，当地时间晚上 8:00(GMT+2)斯德哥尔摩的一个朋友的在线状态与当地时间晚上 8:00(GMT-4)纽约州的另一个朋友的在线状态进行比较。
- 因为我的大多数朋友都来自巴基斯坦，默认时区被假定为 GMT+5，对于我碰巧记得的朋友，我相应地调整了他们的时区。然而，由于时区设置错误，汇总数据中可能存在一些不准确性(约 5-10%)。
背景信息已经够多了，现在让我们来看看第一周和第二周在线朋友的平均数量。图 1 描绘了一天中不同时间线上朋友的数量。
Figure 1: Average number of online friends over the day for two different weeks
- 这两周的曲线高度略有不同，但这两周曲线模式在一天不同时间的变化非常一致。
- 两周后曲线的峰值非常一致。
如果您想知道如何创建这些图表，请继续，您可以在本文的下一部分中了解到这一点。
从图 1 中可以明显看出，一周的数据通常是对模式的足够好的估计，我们可以继续绘制一周中不同日子的在线朋友数量图。
Figure 2: Number of online friends over the day for whole week
让我们列出图 2 中的一些观察结果:
- 一天中在线朋友数量的模式在一周中不同的日子里惊人地一致，高度几乎没有变化。
- 当地时间下午 5:00 为每个好友记录了最多的在线好友数量。
- 凌晨 4:00 记录的在线好友数量最低(可能是因为大多数好友都在睡觉)。
我们还可以绘制一个条形图，显示一周中不同日子的平均在线时间。
Fig 3: Average number of hours spent online for each day of the week
- 朋友们平均每天在脸书呆 3.5 个小时。
- 周六和周一被发现是一周中最忙的一天，同时也是上网时间最多的一天。
如本节开头所述，出于某些原因，每个用户都被视为基于其当地时间在线。但是，如果我们想知道一天中最大限度地扩大帖子覆盖范围和浏览量的最佳时间，这可能会产生误导性的统计数据。
让我们根据我的时区(恰好是中欧时间)调整每个朋友的时区，然后绘制一天中在线朋友的数量图，以确定一天中发布新帖子的最佳时间。
Figure 4: Number of online friends over the day for whole week (relative to my local time)
到目前为止，我们只使用了用户的汇总数据，但是个人朋友的在线活动呢？我们能从中获得什么样的见解？
假设用户倾向于把查看脸书作为晚上的最后一件事，也是早上的第一件事，我们可以发现在脸书上不活跃的一些非常一致的模式，这表明用户正在睡觉。
让我们画出我的一个朋友的睡眠模式。从图 5 中可以观察到，我的朋友从午夜 12:00 到早上 9:00 平均每天睡眠 6.5 小时，非常规律。
Figure 5: Sleep times of a particular user over days
仅通过可视化一天中用户在线/离线的时间，就可以很好地估计用户的睡眠模式。
很有趣，不是吗？嗯，这是 迷人的 同时也是 令人担忧的 ，让我们意识到这么小的一条信息能透露多少关于我们日常生活的信息。
这些都是通过分析脸书朋友的在线活动得出的见解。如果您有兴趣了解如何使用 Python 创建本节中呈现的图表，可以继续阅读下一节。否则，你可能会有兴趣阅读**【即将推出】分析社交网络图。**
使用 Python 可视化脸书朋友的在线活动
先决条件:
- Python 3.5
- 钟摆
- Matplotlib
- Numpy
为了简单和隐私起见，我不会在这里展示我朋友的真实数据。在本教程中，我们将使用具有相同结构的样本虚拟数据。
让我们看看数据集的结构。收集的数据集是一个简单的字典列表，其中包含“在线朋友列表”和“时间戳”(当列表被获取时)，作为 JSON 存储在一个文件中。**
*[
{
"timestamp": "2019-05-14T08:12:52.313921+02:00",
"friends": [
"John Doe",
"John Roe"
]
},
{
"timestamp": "2019-05-14T08:13:49.449888+02:00",
"friends": [
"Jane Roe",
"Jane Doe"
]
},
...
...
]*
加载数据
让我们加载样本数据并执行一些初始处理:
上面的代码将为我们提供“时间戳”，这是一个(friend_name，times)对的 Python 字典，其中 times 是一个列表，其中的pendulum . datetime对象表示特定朋友在线时记录的时间戳。
下一步是将 dict 转换成嵌套的 dict，每个朋友的时间戳按天分组。下面的代码可以做到这一点:
上述代码执行以下操作:
- 它按照唯一的日期对每个朋友的 pendulum.datetime 对象列表进行分组，并将其转换为(day，times)对的字典，其中 times 是按时间顺序排列的时间戳列表。
- 创建一个名为“activities”的嵌套有序字典，其中包含一个嵌套字典对(day，times)。
结果是以下格式的嵌套字典:
*{
"Jane Roe": {
"We 15/05/19": [DateTime,...,],
"Th 16/05/19": [DateTime,...,]
},
"Jane Doe": {
"We 15/05/19": [DateTime,...,],
"Th 16/05/19": [DateTime,...,]
},
...
}*
给定上面的字典，我们想要做的是创建一个二进制地图来表示每个朋友每天的在线状态。就这么办吧。
如果您觉得它看起来吓人，请不要担心，我还在代码片段后提供了代码摘要。
上面的代码实质上是将每个朋友每天的时间戳列表转换成一个二进制映射，1 表示在线，0 表示离线，时间间隔为 5 分钟，时间间隔为一天 24 小时。它执行以下操作:
- 一天创建 289 个垃圾箱。每个时间戳被放入适当的容器中(但是不计算频率)。
- 创建一个名为 online_maps 的 shape (4，2，289)的 Numpy 数组，因为有 4 个朋友和 2 天。创建阵列的日期按时间顺序排列。
为了熟悉箱的概念，假设我们将一天分成 24 个部分/间隔，即箱，每小时一个部分。如果用户在从午夜 12:00 到凌晨 1:00 的任何时间在线，我们将把第一个 bin 的值设置为 1，否则为 0。对所有 24 个时间间隔执行此操作将产生一个整数列表(我们称之为二进制映射，因为每个整数的值为 0 或 1 ),表示一天中每个小时用户在线或离线的状态。
我们现在有一个非常好的数据格式，即每个朋友几天的活动地图。让我们使用格式化的数据绘制一些图表。
绘图活动状态
我们可以画出一个朋友在不同日子不同时间的在线活动。下面的代码描绘了数据中第一个朋友的活动。
这将绘制下图，其中蓝色竖线表示在线，不在线表示离线。
Figure 6: Online status of a friend across different days of the week
有趣到目前为止？
让我们绘制另一个图表，在这个图表中，我们可以看到一天中不同时间的在线好友总数。
在线好友的数量
通过计算 online_maps 数组第一维的总和，我们可以得到一天中在线的朋友数量，即朋友数量:
*>>> online_count = np.sum(online_maps, axis=0)
>>> print(online_count.shape)
(2, 289)*
下面的代码绘制了一个折线图，显示不同日子在线的朋友数量。为了视觉效果，我们也将高斯平滑应用于每条曲线。
运行上面的代码会产生下面的图:
Figure 7: No. of online friends of different hours of multiple days in the week
利用这组简单的信息，我们还可以做更多的事情(其中一些我将在本系列的后面部分讨论)，但是本教程应该可以帮助您从基础开始。本帖完整代码的 Jupyter 笔记本可以在这里找到。
我希望这篇文章对你的工作和研究有所帮助。如果你有任何问题，欢迎在下面写评论。
阅读相关故事:
- 抓取你脸书朋友的数据
- 分析在线活动和睡眠模式
- ****【即将推出】:分析我的社交网络图
分析招聘数据
无需提及统计学或机器学习就能获得洞察力
Photo by Clem Onojeghuo on Unsplash
上周，我发现自己在和一位同事谈论我们多年来参与的各种分析项目。不用说，对话冒险进入了复杂的统计分析、技术术语和超级有趣的发现的兔子洞。在谈话结束时，我肯定会提起我几年前进行的一项具体分析，该分析揭示了一系列可操作的见解，但没有提到统计意义。没有提到相关或混淆矩阵。没有提到张量流或递归神经网络。
搭建舞台
大约两年前，我的一个朋友参加了一个 SuccessFactors 会议，会上我展示了 SuccessFactors 平台自带的分析工具的强大功能。会议结束后，他吸引了我几分钟的注意力，问我是否可以帮助他解决他遇到的招聘问题。尽管尽了最大努力，他还是很难将高潜力的候选人转化为雇员。这个问题是短期公司战略的一个重大障碍，因为随着公司准备进入一个新产品市场，它需要多样化的技能组合。我已经让他把他的数据发过来，我会看一看。
分析/结果
收到数据后，我决定将招聘过程分成三个独立的部分；申请、申请人和征用。
应用程序
首先，我通过得出“申请完成率”、“申请退出率”和“完成申请的平均时间”等指标来检查申请流程。令我惊讶的是，完成率超过 80%，脱落率极低，提交申请平均需要 3 分钟。问题当然不在于申请过程。
申请人
接下来，我把注意力转向了申请人。哪些来源推动了申请人群体？哪些来源产生了高质量的申请人？拥有庞大的申请人才库大大增加了聘用高质量候选人的可能性。分析的结果产生了混合的结果，因为申请人才库是足够的，但高度倾斜。17 个候选来源中的 3 个驱动了 80%的申请人和大多数优质申请人(即进入招聘经理面试阶段的候选人)。这些结果并不一定会给我指出任何可能解决潜在招聘问题的方向，但我发现了潜在的效率和成本节约因素。这些较小的见解在数据挖掘实践中非常常见。
征用
接下来，我开始更仔细地分析征用数据。首先，我检查了完成一个申请平均需要多长时间(完成时间定义为申请完成日期—申请发布日期)。我特别关注技术需求，因为它们对业务战略有最直接的影响。令我惊讶的是，填写一份申请表平均需要 50 多天。需求与头寸的 1:1 比率)这一指标在过去六个月中没有逐月波动。更令人不安的是，同一时期的招聘率低于 5%。这些结果肯定证实了我朋友提出的担忧，但是是什么导致了这些问题呢？在需求技术候选人往往可以有两个或三个报价，时间是至关重要的。我发现这家公司雇佣了 8 名招聘人员，他们管理着 15 到 25 个职位需求，这不一定是很大的工作量，但是当时间非常重要时，这肯定是问题的一部分。
我把注意力集中在招聘过程本身，以更好地理解是什么导致了如此激烈的时间填补时代。通过候选人漏斗分析，我发现了一个 12 步招聘流程！在候选人进入“评估”阶段之前，阶段之间的候选人退出率和阶段中的平均时间都在正常范围内。在进入这一阶段的候选人中，只有 4%的人进入了下一阶段“创建录用通知”。此外，绝大多数被处分的候选人都自愿终止了申请过程！审查阶段的平均时间从候选人被评估到进入“创建录用通知”阶段平均需要 16 天，实际收到录用通知又需要 6 天。我发现了影响招聘对话率的最后一块拼图。
推荐
最终的罪魁祸首是候选人的满意度。招聘人员不得不处理太多的申请，这意味着他们对候选人的数量负担过重。雇佣更多精通技术的招聘人员可以分散工作量，让候选人更快地接受面试，最终让他们更高效地通过面试。此外，我建议改进申请流程，以减少候选人在获得聘用之前必须经历的阶段。这将加快招聘过程，而不会显著影响候选人的质量。为招聘人员提供培训和标准化的行为评估将有助于消除选拔过程中的偏见，产生质量相似的候选人，并加快候选人的选拔过程。接下来，我建议详细审查聘用创建工作流程，以简化流程并减少候选人收到聘用的天数。最后，出于节约成本的考虑，我建议公司进一步检查来源预算分配，因为绝大多数候选人都来自少数来源。最有效的来源渠道之一，推荐，并没有得到充分利用。
拿走
我们经常偶然发现一些鼓舞人心的文章，这些文章讲述了一些公司利用先进的技术和统计知识来获得对组织产生巨大影响的见解。其中一些分析所需的知识、技能和大量数据简直令人难以置信。这并不意味着不能从单一数据源和使用 Microsoft excel 的精明分析师那里收集可操作的见解。对数据、度量设计、技术诀窍和问题的背景知识有敏锐的理解会大有帮助。
Youtube 数据上的文本分类技术分析
文本分类是自然语言处理所要解决的一个经典问题，它指的是分析原始文本的内容并确定其所属的类别。这类似于有人读了罗宾·夏尔马的书，把它归类为‘垃圾’。它有广泛的应用，如情感分析，主题标记，垃圾邮件检测和意图检测。
今天，我们将进行一项相当简单的任务，使用不同的技术(朴素贝叶斯、支持向量机、Adaboost 和 LSTM)并分析其性能，根据视频的标题和描述将视频分类为不同的类别。这些类别被选择为(但不限于):
- 旅游博客
- 科学与技术
- 食物
- 制造业
- 历史
- 艺术和音乐
事不宜迟，就像一个刚刚开始园艺的中年父亲会说，“让我们把手弄脏吧！”。
收集数据
当处理像这样的定制机器学习问题时，我发现收集自己的数据非常有用，如果不是简单地令人满意的话。对于这个问题，我需要一些关于属于不同类别的视频的元数据。如果你有点笨，我欢迎你手动收集数据并构建数据集。然而我不是，所以我会用 Youtube API v3 。它是由谷歌自己创建的，通过一段专门为我们这样的程序员编写的代码与 Youtube 进行交互。前往谷歌开发者控制台，创建一个示例项目并开始。我选择这样做的原因是，我需要收集成千上万的样本，我发现使用任何其他技术都不可能做到这一点。
注意:Youtube API，像谷歌提供的任何其他 API 一样，基于配额系统工作。根据你的计划，每封邮件每天/每月都有固定的配额。在我的免费计划中，我只能向 Youtube 发出大约 2000 次请求，这造成了一点问题，但我通过使用多个电子邮件帐户克服了这个问题。
API 的文档非常简单，在使用了超过 8 个电子邮件帐户来补偿所需的配额后，我收集了以下数据并将其存储在一个. csv 文件中。如果你想在你的项目中使用这个数据集，你可以在这里 下载 。
Collected Raw Data
注意:你可以自由探索一种被称为网络抓取的技术，它被用来从网站中提取数据。Python 有一个名为 BeautifulSoup 的漂亮库用于同样的目的。然而，我发现在从 Youtube 搜索结果中抓取数据的情况下，对于一个搜索查询，它只返回 25 个结果。这对我来说是一个障碍，因为我需要大量的样本来创建一个准确的模型，而这并不能解决问题。
数据清理和预处理
我的数据预处理过程的第一步是处理丢失的数据。因为丢失的值应该是文本数据，所以没有办法估算它们，因此唯一的选择是删除它们。幸运的是，在总共 9999 个样本中，只存在 334 个缺失值，因此它不会影响训练期间的模型性能。
“视频 Id”列对于我们的预测分析并不真正有用，因此它不会被选为最终训练集的一部分，所以我们没有任何预处理步骤。
这里有 2 个重要的栏目，即— 标题和描述，但它们都是未经加工的原文。因此，为了过滤掉噪声，我们将遵循一种非常常见的方法来清理这两列的文本。这种方法分为以下几个步骤:
- **转换成小写:**执行这一步是因为大写对单词的语义重要性没有影响。“travel”和“Travel”应该被同等对待。
- **删除数值和标点:**数值和标点中使用的特殊字符($，！等等。)无助于确定正确的类别
- **删除多余的空格:**这样每个单词都由一个空格分隔，否则在标记化过程中可能会出现问题
- **记号化成单词:**这指的是将一个文本串分割成一个“记号”列表，其中每个记号是一个单词。例如，句子“我有巨大的二头肌”将转换为[“我”、“有”、“巨大”、“二头肌”]。
- 去除非字母单词和“停用词”:“停用词”指的是 and、the、is 等单词，这些单词在学习如何构造句子时很重要，但对我们的预测分析没有用处。
- **词汇化:**词汇化是一种非常棒的技术，可以将相似的单词转换成它们的基本含义。例如，单词“flying”和“flyed”都将被转换成它们最简单的意思“fly”。
Dataset after text cleaning
“现在文字干净了，万岁！让我们开一瓶香槟庆祝吧！”
不，还没有。即使今天的计算机可以解决世界问题，玩超逼真的视频游戏，但它们仍然是不懂我们语言的机器。因此，我们无法将文本数据原样提供给我们的机器学习模型，无论它有多干净。因此，我们需要将它们转换成基于数字的特征，这样计算机就可以构建一个数学模型作为解决方案。这就构成了数据预处理步骤
Category column after LabelEncoding
由于输出变量(’ Category ')也是分类的，我们需要将每个类编码为一个数字。这叫做标签编码。
最后，让我们关注一下每个样本的主要信息—原始文本数据。为了从文本中提取数据作为特征并以数字格式表示它们，一种非常常见的方法是将它们矢量化。Scikit-learn 库包含用于此目的的“TF-IDF 矢量器”。 TF-IDF (词频-逆文档频率)计算每个词在多个文档内部和跨文档的频率，以识别每个词的重要性。
数据分析和特征探索
作为一个额外的步骤，我决定显示类的分布，以便检查样本数量的不平衡。
此外，我想检查使用 TF-IDF 矢量化提取的特征是否有意义，因此我决定使用标题和描述特征为每个类找到最相关的单字和双字。
# USING TITLE FEATURES# 'art and music':
Most correlated unigrams:
------------------------------
. paint
. official
. music
. art
. theatre
Most correlated bigrams:
------------------------------
. capitol theatre
. musical theatre
. work theatre
. official music
. music video# 'food':
Most correlated unigrams:
------------------------------
. foods
. eat
. snack
. cook
. food
Most correlated bigrams:
------------------------------
. healthy snack
. snack amp
. taste test
. kid try
. street food# 'history':
Most correlated unigrams:
------------------------------
. discoveries
. archaeological
. archaeology
. history
. anthropology
Most correlated bigrams:
------------------------------
. history channel
. rap battle
. epic rap
. battle history
. archaeological discoveries# 'manufacturing':
Most correlated unigrams:
------------------------------
. business
. printer
. process
. print
. manufacture
Most correlated bigrams:
------------------------------
. manufacture plant
. lean manufacture
. additive manufacture
. manufacture business
. manufacture process# 'science and technology':
Most correlated unigrams:
------------------------------
. compute
. computers
. science
. computer
. technology
Most correlated bigrams:
------------------------------
. science amp
. amp technology
. primitive technology
. computer science
. science technology# 'travel':
Most correlated unigrams:
------------------------------
. blogger
. vlog
. travellers
. blog
. travel
Most correlated bigrams:
------------------------------
. viewfinder travel
. travel blogger
. tip travel
. travel vlog
. travel blog# USING DESCRIPTION FEATURES# 'art and music':
Most correlated unigrams:
------------------------------
. official
. paint
. music
. art
. theatre
Most correlated bigrams:
------------------------------
. capitol theatre
. click listen
. production connexion
. official music
. music video# 'food':
Most correlated unigrams:
------------------------------
. foods
. eat
. snack
. cook
. food
Most correlated bigrams:
------------------------------
. special offer
. hiho special
. come play
. sponsor series
. street food# 'history':
Most correlated unigrams:
------------------------------
. discoveries
. archaeological
. history
. archaeology
. anthropology
Most correlated bigrams:
------------------------------
. episode epic
. epic rap
. battle history
. rap battle
. archaeological discoveries# 'manufacturing':
Most correlated unigrams:
------------------------------
. factory
. printer
. process
. print
. manufacture
Most correlated bigrams:
------------------------------
. process make
. lean manufacture
. additive manufacture
. manufacture business
. manufacture process# 'science and technology':
Most correlated unigrams:
------------------------------
. quantum
. computers
. science
. computer
. technology
Most correlated bigrams:
------------------------------
. quantum computers
. primitive technology
. quantum compute
. computer science
. science technology# 'travel':
Most correlated unigrams:
------------------------------
. vlog
. travellers
. trip
. blog
. travel
Most correlated bigrams:
------------------------------
. tip travel
. start travel
. expedia viewfinder
. travel blogger
. travel blog
建模和培训
我们将分析的四个模型是:
- 朴素贝叶斯分类器
- 支持向量机
- Adaboost 分类器
- LSTM
数据集被拆分为训练集和测试集，拆分比例为 8:2。标题和描述的特征被独立计算，然后被连接以构建最终的特征矩阵。这用于训练分类器(除了 LSTM)。
对于使用 LSTM，数据预处理步骤与前面讨论的非常不同。这是它的流程:
- 将每个样本的标题和描述组合成一个句子
- **将组合句子标记为填充序列:**将每个句子转换为标记列表，为每个标记分配一个数字 id，然后通过填充较短的序列并截断较长的序列，使每个序列具有相同的长度。
- 对“Category”变量进行一次性编码
LSTM 的学习曲线如下所示:
LSTM Loss Curve
LSTM Accuracy Curve
分析性能
以下是所有不同分类器的精度-召回曲线。要获得额外的指标，请查看 完整代码 。
在我们的项目中观察到的每个分类器的等级如下:
LSTM > SVM >朴素贝叶斯> AdaBoost
LSTMs 在自然语言处理的多个任务中表现出了优异的性能，包括这个任务。LSTMs 中多个“门”的存在允许它们学习序列中的长期依赖性。深度学习 10 分！
支持向量机是高度鲁棒的分类器，它们尽最大努力寻找我们提取的特征之间的相互作用，但是学习到的相互作用与 LSTMs 不一致。另一方面，朴素贝叶斯分类器认为特征是独立的，因此它的性能比支持向量机稍差，因为它没有考虑不同特征之间的任何交互作用。
AdaBoost 分类器对超参数的选择非常敏感，因为我使用了默认模型，所以它没有最佳参数，这可能是性能差的原因
我希望这对你和对我一样有益。完整的代码可以在我的 Github 上找到。
拥有 3 年以上项目工作经验的初级数据科学家和软件开发人员。高度熟练的机器…
github.com](https://github.com/agrawal-rohit)
再见
使用熊猫分析疾病控制中心(CDC)的癌症数据，第 1 部分
Photo by Miguel Á. Padriñán on Pexels
虽然癌症药物开发和治疗的研究多年来取得了进展，但癌症每年仍夺去成千上万人的生命。尽管如此，随着对数据的访问、计算能力和最先进的机器学习工具的增加，癌症研究进展的潜力继续增长。
在本帖中，我们将探索疾病控制中心癌症数据集。该数据集包括脑瘤、各州癌症类型、种族、年龄等信息。在这篇文章中，我们将探索“大脑旁站点”。’ TXT '数据。
我们从导入熊猫图书馆和阅读开始。TXT '文件转换成熊猫数据帧。每一列都用“|”分隔，所以我们也适当地设置分隔参数“sep”。我们还可以指定需要分析的列，并显示前五行数据，以了解列类型和值:
df = pd.read_csv("BRAINBYSITE.TXT", sep="|")
df = df[['AGE', 'BEHAVIOR', 'COUNT', 'POPULATION', 'SEX', 'YEAR', 'SITE"]]
print(df.head())
“计数”列包含一些我们可以删除的缺失值:
df = df[df['COUNT'] != '~']
df.reset_index(inplace=True)
print(df.head())
为了开始我们的分析，我们可以生成“计数”列的直方图，以可视化所有类别的肿瘤分布:
import seaborn as sns
import matplotlib.pyplot as plt#settings for the histogram plot
sns.set(font_scale = 2)
plt.ylim(0, 80)
plt.xlim(0, 10000)df['COUNT'] = df['COUNT'].astype(int)
df['COUNT'].hist(bins=1000)
Distribution in Tumor Counts
我们还可以在重叠图上查看女性和男性肿瘤数量的直方图:
#define female and male data frames
df_female = df[df['SEX'] == 'Female']
df_male = df[df['SEX'] == 'Male']#overlay histograms
df_female['COUNT'].hist(bins=1000)
df_male['COUNT'].hist(bins=1000).set_title('Male and Female')
Overlay of distributions in tumor counts for males and females
我们可以进行更细致的分析，观察给定年份(假设“年份”=2004 年)女性和男性的肿瘤数量分布:
df_female = df[df['SEX'] == 'Female']
df_female = df_female[df_female['YEAR'] == '2012-2016']
df_female.loc[:, 'COUNT'] = df_female['COUNT'].astype(int)
df_female['COUNT'].hist(bins=1000).set_title('Female ')
Distribution in tumor counts for females between 2012–2016
df_male = df[df['SEX'] == 'Male']
df_male = df_male[df_male['YEAR'] == '2012-2016']
df_male.loc[:, 'COUNT'] = df_male['COUNT'].astype(int)
df_male['COUNT'].hist(bins=1000).set_title('Male')
Distribution in tumor counts for males between 2012–2016
接下来，我们可以生成“原油价格”与“年份”的散点图。在这里，我们还将调整散点图点的大小，以便大小和色调与“计数”成比例。我们还过滤数据帧以获得对应于 0.1 以上的“部位”=“间变性星形细胞瘤”、“行为”=“恶性”和“粗制率”的行。对于女性，我们有:
df_female = df[df['SEX'] == 'Female']
df_female = df_female[df_female['SITE'] == 'Anaplastic astrocytoma']
df_female = df_female[df_female['BEHAVIOR'] == 'Malignant']#remove baseline value present in each year
df_female = df_female[df_female['CRUDE_RATE'] > 0.1]sns.scatterplot(df_female["YEAR"], df_female["CRUDE_RATE"], sizes = (1000, 1500), size = df_female["COUNT"], alpha = 0.8,hue = df_female["COUNT"])
Scatter plot of ‘CRUDE_RATE’ vs ‘YEAR’ of Anaplastic astrocytoma in females
我们也可以看看这个雄性的图:
df_male = df[df['SEX'] == 'Male']
df_male = df_male[df_male['SITE'] == 'Anaplastic astrocytoma']
df_male = df_male[df_male['BEHAVIOR'] == 'Malignant']#remove baseline value present in each year
df_male = df_male[df_male['CRUDE_RATE'] > 0.1]
sns.scatterplot(df_male["YEAR"], df_male["CRUDE_RATE"], sizes = (1000, 1500), size = df_male["COUNT"], alpha = 0.8,hue = df_male["COUNT"])
Scatter plot of ‘CRUDE_RATE’ vs ‘YEAR’ of Anaplastic astrocytoma in males
如你所见，2004 年至 2016 年间，女性和男性间变性星形细胞瘤的粗发病率都有所上升。还值得注意的是，男性患这种脑瘤的比率高于女性。尝试生成其他肿瘤类型的散点图，看看女性和男性之间是否有任何有趣的差异。
接下来，我们可以从这些数据列中生成一些统计数据。我们可以定义一个函数，输出分类变量的唯一值集，并计算该值在数据中出现的次数:
from collections import Counter
def get_unqiue_values(feature):
print("{} Unique Set: ".format(feature), set(df[feature]))
print("{} Count: ".format(feature), dict(Counter(df[feature])))get_unqiue_values('SEX')
当我们使用“性别”字段调用“get_unqiue_values”函数时，输出显示每个“性别”类别的唯一值列表以及每个值在数据中出现的次数。
Uniques set of ‘SEX’ values and their frequencies
我们可以对“行为”做同样的事情:
get_unqiue_values('BEHAVIOR')
输出是“行为”类别的唯一值的列表，以及每个值在数据中出现的次数。
Uniques set of ‘BEHAVIOR’ values and their frequencies
对于“站点”类别:
get_unqiue_values('SITE')
输出显示了“SITE”列的唯一值列表以及每个值在数据中出现的次数。
Unique set of ‘SITE’ values and their frequencies
与该数据集中的其他分类变量相比，存在相对更多的唯一“位置”值。我们可以通过查看五个最常见(最频繁出现)的值来缩小范围:
from collections import Counter
def get_unqiue_values(feature):
print("{} Unique Set: ".format(feature), set(df[feature]))
print("{} Count: ".format(feature), dict(Counter(df[feature]).most_common(5)))
get_unqiue_values('SITE')
Most common ‘SITE’ values
将这些数据可视化可能比盯着字典更有用一些。为此，我们可以更新函数“get_unqiue_values ”,使其返回 10 种最常见脑瘤的字典:
def get_unqiue_values(feature):
print("{} Unique Set: ".format(feature), set(df[feature]))
print("{} Count: ".format(feature), dict(Counter(df[feature]).most_common(5)))
result = dict(Counter(df[feature]).most_common(10))
return result
接下来，我们迭代“站点”频率字典，并将键和值存储在数据框中:
key_list = []
value_list = []
for key, value in get_unqiue_values('SITE').items():
key_list.append(key)
value_list.append(value)
site_df = pd.DataFrame({'SITE': key_list, 'Count':value_list} )
最后，我们使用 seabon 定义一个条形图对象，并设置 x 轴标签:
ax = sns.barplot(x=site_df.SITE, y=site_df.Count)
ax.set_xticklabels(site_df.SITE, rotation=30)
The frequency of tumor types
在‘BRAINBYSITE’有更多值得探索的地方。TXT '数据，但现在我们将结束我们的分析。在接下来的几篇文章中，我将继续探索疾病控制中心癌症数据集中提供的一些其他数据集。在那之前，你可以自己在其他数据集上重复这个分析。这篇文章中显示的代码可以在 GitHub 上找到。祝好运，机器学习快乐！
分析美国地震以确定可能的地震风险。
Photo by Sunyu Kim on Unsplash
曾经想过，如果有一个预测系统，可以在地震袭击我们脚下的地面之前准确地预测地震，死亡人数和几座摩天大楼、住宅、酒店的破坏将大幅减少，因此我们可以逃离自然的最大灾难之一，地震将不再成为一个人死亡的原因。
假设，你去了美国太平洋沿岸的一个海滨城市，并入住了一家酒店。现在，你正在欣赏海景，突然你感觉到一个脆弱的震动，酒店，可能是一座摩天大楼像玩具一样倒塌了，现在你想知道如果你知道在这些地震严重的地方哪些酒店需要仔细检查和升级抗震设备，你不会只受用户评级的影响而选择这家酒店入住。
但是，通过当今非常强大的数据科学工具的进步，有可能确定发生地震的可能性很高的地方，通过那个地方，人们可以很容易地获取附近商业场所的数据，包括商业场所、公寓或任何东西，以测量风险因素和结构强度，为它们配备抗震设备并升级它们，从而确保居民和游客的安全。
那么，我们要做什么？高层次概述:
所以，现在，问题对我们来说是清楚的。我们需要分析地震记录，必须确定地震发生率较高的地点，并确定因其位置而可能处于危险中的酒店。考虑到这个问题，让我们采取有条不紊的方法来解决这个问题。我们将按照约翰·罗林斯博士提出的方法来解决这个问题。
首先，让我们从某个数据提供商那里获取一个包含美国地震数据的数据集，开始处理这些数据。在获取阶段完成后，让我们继续清理获取的数据集，方法是删除不必要的列，用有意义的内容替换空值或无效值，并通过用不规则范围规范化值来使所有列同等重要。
现在我们分析数据框的不同列以及它们与其他列的关系。简而言之，为了便于分析，我们现在试图找到列值之间的关系。在这个过程中，我们使用了几个 Python 可视化库和包来支持我们的发现。
现在，我们对整个数据框进行分组，以搜索出现频率独特的位置。
接下来，我们将继续收集记录地震事件最频繁的地方的地理编码纬度和经度。我们将把这个纬度和经度数据传递给我们的街道级信息提供程序，以获取一个包含不同商业场所信息的结果集，包括记录的场所附近的酒店信息。
选择地震和街道级数据提供商并请求数据:
在我们深入分析之前，我们需要大量的数据和信息。在数据科学领域，选择可靠且价格合理的数据提供商是至关重要的。我们必须确保提供给我们的数据是高效、可靠、使用灵活、可修改和定制的。事实上，如果是人工审核，而不是自动生成，准确性和可靠性会提高。
所以，考虑到所有的事情，我决定和美国地质调查局(USGS)一起收集地震数据。 USGS 提供了非常直观、易用、可靠的 API 和 web 门户服务，在输出格式、指定感兴趣的区域等方面提供了灵活性。美国地质勘探局是一个政府运营的研究中心，他们提供的数据是免费的，非常可靠，因为大多数数据在注册前都经过了人工审查。API 请求链接不需要任何身份验证。
因此，我们指定我们感兴趣的地区为美国，选择我们感兴趣的时间段为 2010 年至 2019 年，输出格式为逗号分隔值(CSV ),然后点击搜索以获得结果。我们可以使用 API 来请求数据。请求 URL 的格式如下所示。
import requests
usgs_request_url = ‘[https://earthquake.usgs.gov/fdsnws'](https://earthquake.usgs.gov/fdsnws')
usgs_request_url+=’/event/1/query.geojson?starttime=2010–04–20%2000:00:00&endtime=2019–04–27%2023:59:59'
usgs_request_url+=’&maxlatitude=50&minlatitude=24.6&maxlongitude=-65&minlongitude=-125&minmagnitude=2.5&orderby=time&limit=2000'requests.get(usgs_request_url).json()
requests
我已经下载了一个格式化的 CSV 数据框架来使用。可以按照链接:https://S3 . Amazon AWS . com/swastiknath group/seismic _ results . CSV下载文件，在自己的笔记本或其他任何 python 环境下使用。
获取数据后，我们使用 pd.read_csv()函数从 CSV 文件创建一个 Pandas Dataframe，并使用 pd.drop()函数删除所有不必要的列，该函数提供要作为 list 对象删除的列，作为 drop 函数的输入参数。
我还决定使用 Foursquare API 来获取美国各地商业场所的街道信息，包括我感兴趣的地方。他们的 API 使用起来非常便宜，有很好的文档记录并且可靠。每次连接到 Foursquare 服务器并请求数据时，您都必须验证您的 API 请求 URL。身份验证的格式如下:
#credentials for communicating with Foursquare API:
CLIENT_ID= 'your client id'
CLIENT_SECRET = 'your client secret'
VERSION = integer date of version.
一旦我们提取完地震数据，我们将定义一个 python 函数，将 USGS 地震数据集返回的纬度和经度列表作为输入，并使用该数据搜索易受地震影响的地方附近的酒店。自定义函数 whats _ nearby 也从 JSON 响应和 Pandas 数据框中检索数据，以更好地表示这些值。
#custom function for requesting street-level business information #from foursquarefrom progressbar import ProgressBar, Bar, Percentage
pbar = ProgressBar(widgets=[Percentage(), Bar()], maxval=300).start()
import requests
def whats_nearby(latitudes, longitudes, radius = 15000):
nearby_places_list1 = []
status = 'EXECUTION STARTED, PLEASE WAIT ->'
for lat, lng in zip(latitudes, longitudes):
for i in range(len(latitudes)):
pbar.update(i+1)
pbar.finish()
query_url = '[https://api.foursquare.com/v2/venues/search?&query=hotel&client_id={}&client_secret={}&v={}&ll={},{}&radius={}&limit={}'.format(](https://api.foursquare.com/v2/venues/search?&query=hotel&client_id={}&client_secret={}&v={}&ll={},{}&radius={}&limit={}'.format()
CLIENT_ID,
CLIENT_SECRET,
VERSION,
lat, lng,
15000,
100)
response = requests.get(query_url).json()['response']['venues']
nearby_places_list1.append([(r['name'],
r['location']['lat'],
r['location']['lng'],
r['location']['formattedAddress']) for r in response])
vulnerable_places = pd.DataFrame(item for nearby in nearby_places_list1 for item in nearby)
vulnerable_places.columns = ['Venue', 'Venue Latitude', 'Venue Longitude', 'Address']
return vulnerable_places
分析变量和它们之间的关系曲线:
现在让我们尝试理解特征集的变量之间的关系，以理解单个变量对地震测量的影响。
想象一下，尽管机器学习、监督、非监督和深度强化学习过程取得了很大进步，但是什么一直在机器学习模型和完全准确的地震预测之间建立障碍？自然和由自然生成的数据是非常不规则的，并且不能被精确地映射成数学关系，尽管自然的每个元素都是由算法驱动的，但是它们理解起来非常非常复杂，并且对于机器学习模型来说是不规则的、不可靠的。每个训练有素的机器学习模型都试图以最简单的方式将特征集与目标细胞进行映射。我们建立了一个深度神经网络来正确预测地震，但仍然没有预期的那么好。
让我们想象一下地震能量波的均方根速度与震级之间的关系。我们使用 Seaborn 包，它是基于 Matplotlib 构建的，可以自动生成具有最高可信度的线性回归图。这是一个很大的特点，它根据回归线显示了曲线的偏差。在均方根与幅度的关系中，我们看到两者之间存在线性关系，尽管线的两侧存在极小的偏差。由于自然的不规则性，这些偏差可以忽略。
sns.regplot(x= earth_quake_encoded['rms'], y= earth_quake_encoded['mag'], data = earth_quake_encoded, color = 'purple')
Regression plot between Root Mean Square Speed and Magnitude
现在，让我们想象地震波破裂的深度与震级之间的关系。由于自然的不规则性导致的偏差量很小的回归线。
import seaborn as sns
reg_plot_depth_mag = pd.DataFrame(columns =['Depth', 'Magnitude'])
reg_plot_depth_mag['Depth'] = earth_quake_encoded['depth'].values
reg_plot_depth_mag['Magnitude'] = earth_quake_encoded['mag'].valuessns.regplot(x=reg_plot_depth_mag['Depth'], y = reg_plot_depth_mag['Magnitude'], data = reg_plot_depth_mag, color = 'blue')
Regression plot between Magnitude and Depth.
现在是时候将从记录到震中的水平距离与震级之间的关系曲线可视化了。回归图显示了具有高偏差的线性回归线，基本上意味着两个变量之间的非线性。
sns.regplot(x= earth_quake_encoded['dmin'], y = earth_quake_encoded['mag'], data = earth_quake_encoded, color = 'crimson')
Regression plot between Horizontal distance vs. Magnitude
现在让我们通过可视化来验证数据集。让我们看一下测量星等误差与星等之间的散点图，以获得对数据集的信心，我们决定也绘制回归线。结果是我们所期望的。回归线就其斜率而言是负的，并且具有最小的偏差，这确实意味着，如果测量量值的误差增加，则实际量值减小。
sns.regplot(x = earth_quake_encoded['magError'], y=earth_quake_encoded['mag'], color = 'orange')
Negative Regression plot between Magnitude and error in measuring magnitude
我们现在想看一下时间框架，基本上是相对于时间连续绘制震级。我们使用 Pandas Series.from_csv 创建一个系列对象，并使用 Matplotlib 的 pyplot 脚本层直接绘制它。在绘制之前，我们确保正确理解了每一列的数据类型。脚本层足够智能，可以相对于选定的时间段以较小的比例(为了更好的可视化)绘制幅度。
time_vs_mag_bar = pd.DataFrame(columns = {'time', 'mag'})
time_vs_mag_bar['time'] = usa_earthquakes['time']
time_vs_mag_bar['mag'] = usa_earthquakes['mag']
time_vs_mag_bar.to_csv('time_vs_mag.csv')dates = mpl.dates.date2num(time_vs_mag_bar['time'].values)
plt.plot_date(dates, usa_earthquakes['mag'])
Time-frame representation of magnitude vs. time
from pandas import Series
series = Series.from_csv('usa_earth_quakes.csv', header=0)
series.plot(style='-')
plt.show()
Pandas Series plotting of continuous time-frame plotting of magnitude(in a smaller scale) vs. time.
我们现在希望在三维空间中可视化数据集，以便更好地反映沿每个轴的变化。因此，我们使用 Matplotlib 的 mpl_toolkits.mplot3d 包的 Axes3D 函数绘制深度、均方根、震级和记录站之间的角度间隙、记录事件所需的站数和震级的三维示意图。
from mpl_toolkits.mplot3d import Axes3D
fig = plt.figure(figsize = (8, 6))
ax = fig.add_subplot(111, projection = '3d')
ax.set_xlabel('Depth')
ax.set_ylabel('RMS')
ax.set_zlabel('Magnitude')
ax.scatter(usa_earthquakes['depth'], usa_earthquakes['rms'], usa_earthquakes['mag'])
plt.show()
3-dimensional representation of Depth vs. RMS vs. Magnitude
fig = plt.figure(figsize = (8, 6))
ax = fig.add_subplot(111, projection = '3d')
ax.scatter(usa_earthquakes['gap'], usa_earthquakes['nst'], usa_earthquakes['mag'])
ax.set_xlabel('Gap')
ax.set_ylabel('Nst')
ax.set_zlabel('Magnitude')
plt.show()
3-dimensional representation of angular gap vs. number of stations vs. magnitude
看看代表美国易受地震影响地区的不同地图:
从检索到的地震数据中，我们获取纬度、经度和地名(使用 Regex 模式获得)。我们使用一个叫做 Folium 的 python 映射库。Folium is 使用 Leaflet.js 创建高度交互式的地图。因此，此时我们在以美国为中心的地图中绘制纬度和经度，并使用一个名为 HeatMap 的 folio 插件来生成一个热图，沿地图显示整个美国的漏洞频率。由于是交互式的，热图会随着缩放级别而变化。我们还可以使用另一个名为 FastMarkerCluster 的功能，以一种整洁的方式可视化地震重灾区。FastMarkerCluster 也是交互式的。它们随着缩放级别折叠和展开。
heat_quake_map = folium.Map(location=[usa_earthquakes['latitude'].mean(),
usa_earthquakes['longitude'].mean()],
zoom_start=4)
latlons = usa_earthquakes[['latitude', 'longitude']].values.tolist()
from folium.plugins import HeatMap
HeatMap(latlons).add_to(heat_quake_map)
heat_quake_map
Heat map showing earthquake heavy regions around USA.
from folium.plugins import MarkerCluster
from folium.plugins import FastMarkerCluster
usa_quake_map = folium.Map(location=[usa_earthquakes['latitude'].mean(),
usa_earthquakes['longitude'].mean()],
zoom_start=4)usa_quake_map.add_child(FastMarkerCluster(usa_earthquakes[['latitude', 'longitude']].values.tolist()))
usa_quake_map
Map showing the clusters of recorded earthquakes and their location in USA.
分析酒店数据:
前面我们看到了如何使用 Foursquare 获取街道级别的数据。因此，在我们获取完酒店数据后，现在是分析从 Foursquare API 返回的附近酒店数据集的时候了。现在，我们删除空值，只需使用’就可以按唯一的酒店名称对值进行分组。数据框的酒店名称列上的“unique()”方法。我们通过简单地使用“获得每个记录位置附近酒店的发生频率(“易受攻击事件”)。“value_counts()”方法。首先，我们返回显示酒店及其易受攻击事件的数据框。
num = fx['Venue'].value_counts()
vulnerables = num.rename_axis('business_name').reset_index(name='Vulnerable Occurrence')
vulnerables
Pandas Dataframe showing the hotel names returned by Foursquare API which are frequently nearby to the recorded earthquake events.
Hotels with their respective Vulnerable occurrences in a horizontal bar chart.
needs_attention = vulnerables
needs_attention.set_index('business_name', inplace = True)
needs_attention.plot(kind = 'barh', width = 0.8, edgecolor = 'black', color = 'tomato')
fig = plt.figure(figsize = (20, 10))
plt.show()from folium.plugins import MarkerCluster
from folium.plugins import FastMarkerCluster
vulnerable_business_map = folium.Map(location = [usa_earthquakes['latitude'].mean(),
usa_earthquakes['longitude'].mean()], tiles = 'OpenStreetMap',
zoom_start=4)
vulnerable_business_map.add_child(FastMarkerCluster(fd_grouped[['Venue Latitude', 'Venue Longitude']].values.tolist()))folium.Circle([usa_earthquakes['latitude'].mean(),usa_earthquakes['longitude'].mean()], radius=1500000, color='red', opacity = 0.6, fill=False).add_to(vulnerable_business_map)
folium.Circle([usa_earthquakes['latitude'].mean(),usa_earthquakes['longitude'].mean()], radius=2000000, color='red', opacity = 0.6, fill=False).add_to(vulnerable_business_map)
folium.Circle([usa_earthquakes['latitude'].mean(),usa_earthquakes['longitude'].mean()], radius=2500000, color='red', opacity = 0.6, fill=False).add_to(vulnerable_business_map)from folium.plugins import HeatMap
latlons = usa_earthquakes[['latitude', 'longitude']].values.tolist()
HeatMap(latlons).add_to(vulnerable_business_map)
vulnerable_business_map
Hotels that requires attention and nearby to the earthquake events on the heat map validating the result.
聚集酒店:
聚类是一种无监督的学习过程，它基于某些相似性分割或聚集未标记的数据集。因此，这里我们将使用 KMeans 作为一种聚类算法，它功能强大但使用简单。它基本上从数据集中选取 k 个质心，并根据距质心的平均 Minkowski 距离将它们动态分配到一个聚类中，直到算法收敛。
这里，我们将数字 5 指定为要将数据分割成的簇的数量。我们像这样使用 ScikitLearn 包:
from sklearn.cluster import KMeans
cluster_list = fx.set_index('Venue')
cluster_list
num_clusters = 5
clusters = cluster_list.drop('Address', 1)
kmeans_cls = KMeans(n_clusters = num_clusters, random_state = 0).fit(clusters)
kmeans_cls.labels_
我们现在删除 hotels 数据集的标签，并将其作为输入提供给 KMeans，以获得聚类标签。最后，我们将聚类标签与它们各自的酒店名称连接起来，并使用颜色索引将它们绘制在叶子地图上，以便更好地可视化聚类。
cluster_list['Cluster Labels'] = kmeans_cls.labels_
cluster_list.reset_index(inplace = True)
cluster_list
Dataframe showing the hotels, their location data and above all, cluster labels.
import matplotlib.cm as cm
import matplotlib.colors as colors
lat_usa = 37.09
lng_usa = -95.71
mapping_earthquake_cluster = folium.Map(location = [lat_usa, lng_usa], zoom_start = 4)
cluster_in = [int(i) for i in cluster_list['Cluster Labels']]
x = np.arange(num_clusters)
ys = [i + x + (i*x)**2 for i in range(num_clusters)]
colors_array = cm.rainbow(np.linspace(0, 1, len(ys)))
rainbow = [colors.rgb2hex(i) for i in colors_array]markers_colors = []
for lat, lon, poi, cluster in zip(cluster_list['Venue Latitude'], cluster_list['Venue Longitude'], cluster_list['Venue'], cluster_in):
label = folium.Popup(str(poi) + ' Cluster ' + str(cluster), parse_html=True)
folium.CircleMarker(
[lat, lon],
radius=10,
popup=label,
color=rainbow[cluster-1],
fill=True,
fill_color=rainbow[cluster-1],
fill_opacity=0.7).add_to(mapping_earthquake_cluster)
mapping_earthquake_cluster
Folium Map showing the clustered hotels.
商业方法:
从商业的角度来说，我们可以说，这个项目将有助于政府城市管理部门，如市政公司或不同的测量师。另一方面，该项目对于不同酒店或企业的利益相关者来说非常有用，以确定他们的酒店或企业是否需要注意抗震设备的升级或结构强度的评估。
要查看完整的代码和它是如何实现的，请点击这里查看我的谷歌合作笔记本。
结束。和平！
“我们是我们的思想造就了我们；所以小心你的想法。文字是次要的。思想活着；他们旅行很远。”
~斯瓦米·维威卡难达。
卫报食谱(第 1 部分):探索性数据
第 1 部分:探索性数据分析——卫报食谱
The Guardian Food Recipes by TheGuardian
介绍
《卫报》每天都会公布名厨的食谱。报纸的美食版为美食爱好者、业余厨师、饥饿的学生或疲惫的父母提供了新菜肴的灵感。它邀请业余厨师探索无限的烹饪世界，接触新的食材和风味，学习烹饪技术或为庆祝盛宴寻找灵感。简而言之,《卫报》的美食栏目吊起了每个读者的胃口。
然而，报纸上的烹饪文章揭示的不仅仅是食物成分。它们影响着饮食文化、消费者购买、读者饮食和传统的家庭食谱。此外，他们还充当年轻厨师的跳板，并维系着值得信赖和喜爱的厨师群体的建立。
我对网络搜集的好奇心以及对数据分析和食物的热情，自然而然地吸引我访问和探索这些报纸文章的内容。《卫报》是我开始研究时想到的第一份报纸，因此让我得以深入了解英国烹饪界。这篇文章的目的是呈现一些从 2009 年 1 月 1 日到 2019 年 6 月 16 日出版的《卫报》食物食谱的探索性数据分析。
如果你觉得有点饿，想在阅读的同时吃点零食，你可以在这里找到所有有趣的食谱！
British Food Recipes — Established Belfast’s poached eggs, Jeremy Fox’s beet gazpacho (photo by Lizzie Mayson), Gill Meller’s roast chicken (photo by Romas Foord), Rebecca Oliver’s smoked mackerel sandwich (photo by Jean Cazals), River Cottage Strawberries (photo by Lizzie Mayson)
数据和方法
在开始我的分析之前，我首先需要收集数据！这是通过使用美汤和请求的网页抓取完成的。因为这是我第一次体验，所以不像烤馅饼那么简单。我花了几个小时通读食谱页面的 HTML 代码，查看每个元素的存储位置。我很幸运地找到了其他媒体成员的巨大帮助，尤其是来自克里·帕克的教程和 Ritvik Kharkar 的文章。
我总共收集了 7，361 篇食谱文章，最早的一篇可以追溯到 2009 年 1 月 3 日，最近的一篇是 2019 年 6 月 16 日。在浏览食谱网址时，我特别注意检索以下七条信息:
- 出版日期
- 文章标题，或者，如果你喜欢，食谱标题
- 文章副标题，一句话描述食谱
- 文章内容、配方文本和配料表
- 作者姓名
- 评论数量和文章被分享的次数
有些文章没有关于上述要素的信息，因此，如果没有找到一个要素，将写下“无”。第二步，我删除了没有作者和字幕的食谱。
分析食谱
以下图表说明了为更好地理解《卫报》的烹饪遗产所采取的分析步骤。在此，我提出以下初步分析:
- 每年每月的食谱数量
- 每年最著名的五位厨师(根据出版的食谱数量)
- 每年最著名的种类(根据食谱数量)
- 为特殊场合发布的食谱的百分比
- 十个最常分享的食谱
- 评论最多的十种食谱
每年每月的食谱数量
下面的条形图显示了 2009 年至 2019 年期间每年发布的食谱数量，以及全年每个月的食谱数量。有趣的是，在 2012 年之后，食品食谱的出版数量有所上升，在 2013 年达到顶峰，共有 882 种食谱。此后，这一数字又略有下降，但从那时起，平均每天发布 2 个食谱。
你对 2019 年的猜测是什么？厨师们会有创意，发表比 2013 年还要多的作品吗？
The bar charts show the number of recipes published each year and each month between 2009 and 2019.
年度五大名厨
看看对《卫报》食谱贡献最大的五位作者，就能洞察厨师们的受欢迎程度。这一年有变化吗？有没有一个厨师得到了更多的关注，并因此获得了在自己的专栏中发表文章的特权？下面的交互图显示了每年发布食谱数量最多的五位厨师。
乍看之下，奈杰尔·斯莱特和约塔姆·奥托林吉是最负盛名的美食食谱作者，这是无可争议的。自 2009 年以来，两位厨师不断为卫报烹饪遗产做出贡献，影响人们购买、准备、烹饪、分享和食用食物的方式。更糟糕的是注意到奈杰尔·斯莱特可能像煎锅上的爆米花一样忙碌，因为他每年贡献超过 100 份食谱。因此，他超越了所有其他厨师！
另一个有趣的点是 2013 年后投稿作者的转变。从 2009 年到 2012 年，厨师丹·莱帕德和休·费恩利-惠汀斯塔定期为《卫报》撰稿。然而，自 2013 年以来，他们的贡献已经停止，让年轻的厨师接手 Top Banana。去年，2018 年，今年，米拉·索达，安娜·琼斯和汤姆·亨特加入了前五名。
The stacked bar chart shows the most popular chefs from 2009 to 2019 (per percentage of the total recipes published each year).
一些快速的事实关于上面的一些厨师！
丹·莱帕德因他惊人的糕点食谱而闻名！如果你正在寻找一个入口即化的食谱，你可能会喜欢丹的食谱。他教授所有烘烤最漂亮、最奢华的糕点、蛋糕和面包的技巧。休·费恩利-惠汀斯托尔因他是河边小屋的跑步者而闻名，也因他发起的与食品生产相关的环境问题运动而闻名。如果你想真正发现用季节性的、当地的和可持续的英国原料准备的典型的英国食物，那么休·费恩利-惠汀斯托尔的食谱将会适合你！
自 2013 年以来，与《卫报》分享美食食谱的厨师名人不断更换。Meera Sodha 因其庆祝印度蔬菜的菜肴而闻名。如果你正在寻找一个去东方目的地的逃脱之旅，Meera 的食谱会让你尝一尝。安娜·琼斯因其蔬菜烹饪而闻名，这种烹饪尊重当地产品的季节性。她透露了在伦敦哪里可以找到最美味的香草作为午餐。最后，汤姆·亨特以尊重自然和人类的多样性而闻名。如果你不知道如何利用你的西兰花茎、你的棕色香蕉或吃剩的面包，那么汤姆的建议会让你大吃一惊。
Highest contributing and new chefs of the Guardian: Nigel Slater, Yotam Ottolenghi, Meera Sodha, Hugh Fearnley-Whittingstall, Anna Jones, Dan Lepard, Tom Hunt (Pictures from their own Instagram account).
年度最佳类别
看看每年的顶级类别，可能会发现与贡献最高的厨师有一些相似之处。这是一些厨师拥有个人食谱专栏的结果，因此他们的大部分食谱都是这样归类的。
在您与下图互动之前，我发现特别有趣的是，有些类别反映了一些全球食品和健康趋势。例如，Yotam Ottolenghi 在 2008 年至 2010 年间写的“新素食主义”类别可能引领了素食主义的当前趋势，以及越来越多的英国人减少饮食中的动物产品。致力于甜食的类别，如丹·勒帕德的“如何烘焙”和鲁比·坦多的“鲁比烘焙”直到 2015 年都在前十名之列。虽然对烘焙的热爱并没有消失，英国烘焙大赛的食谱仍然刊登在《卫报》上，但随着英格兰公共卫生部努力解决肥胖危机，它们的重要性可能在某种程度上有所下降。
The stacked bar chart shows the most popular categories from 2009 to 2019. (Per percentage of the total recipes published each year).
英国的节日和食物
像圣诞节或复活节这样的节日与食物有着内在的联系。这些场合的菜肴都是精心挑选的。一些家庭会提前几天计划菜单，尝试一些新的食谱以及准备传统的家庭食谱！这些场合为厨师、餐馆、食品公司和食品零售商提供了一个激发业余爱好者灵感和诱惑你贪吃的绝佳机会！
下面的图表显示了各种场合的食谱(情人节、复活节、万圣节、篝火之夜、圣诞节、新年)在每年和每月公布的食谱总数中所占的百分比。
The bar chart shows the percentage of recipes for each occasion for each year, each month from 2009 to 2019. (Per percentage of the total recipes published each year, each month)
有趣的是，圣诞节的食谱在 10 月份左右开始出版，在 12 月份占所有食谱的 30%。篝火之夜，也被称为盖伊·福克斯之夜，是英国和一些英联邦国家的传统纪念活动。这个节日在 11 月 5 日庆祝，全国各地都会燃放大型烟花和篝火。通常，家人和朋友会聚在一起享受奶酪土豆、奶油南瓜饼、松脆的太妃糖苹果或粘稠的巧克力棉花糖！
British Festivities Recipes — Liam Charles’ shortbread (photo by Yuki Sugiura), Joe Trivelli’s traditional Easter roast (photo by Jean Cazals), Yotam Ottolenghi’s deviled eggs (photo by Louise Hagger), Yotam Ottolenghi’s alternative Christmas recipes (photo by Louise Hagger), Nigel Slater’s trout tartare, turkey, pumpkin, quince pie (photo by Jonathan Lovekin).
分享最多的 10 道菜，评论最多的 10 道菜
看看十个分享最多的和十个评论最多的食谱，可以让你对读者的兴趣或引起读者争议的文章有一个有趣的了解。下表列出了分享和评论最多的食谱以及它们背后的作者。
毫不奇怪，最著名的厨师也是最受欢迎的。显然，读者喜欢分享他们的圣诞食谱！读者也喜欢学习技巧和诀窍，因此，“如何烹饪”食谱也很受欢迎！另外，请注意“Yotam Ottolenghi 的学生超便宜食谱”的高分享数！学生知道在哪里可以找到最好的食物，对吗？
The table shows the ten most shared recipes (by number of shares)
评论最多的食谱显示了读者对一篇文章的关注程度。这说明文章给了读者思考的食粮！看看下面的食谱和它们的评论，看看食物如何密切影响我们的情绪，有时会引起钦佩、惊讶或愤慨。
The table shows the ten most commented recipes (by number of comments)
结论
该结束了！:-)这篇文章探讨了 2013 年 1 月至 2019 年 6 月发表在《卫报》上的大约 7361 种食物食谱。探索这些食谱并更多地了解著名的英国厨师真是令人愉快。
如果你有兴趣看我的预处理步骤，我是如何挑选，清洗，切碎，调味并最终整理众多数据的，你可以看看我在这里的 GitHub 库！
菜单上的下一道菜是什么？
下一餐(《卫报》食品文章第二部分)将会带来更多你喜欢向其学习的厨师的风味，并展示每个季节你应该准备的食物。我将向您展示我的自然语言处理(NLP)分析的结果。做好准备，看看被切碎的食谱吧！
感谢阅读！奥德
洗碗之前！
我希望你喜欢这篇文章，如果你想留下任何意见，建议或想法，我真的很感激！所以，请随时伸出手来！
也非常感谢 Medium 的成员分享他们的工作，并提供有用的教程，总是让我对有用的工具和方法有新的见解！
PostScriptum:请注意，通过查看《卫报》的美食食谱，我只是为你提供了英国美食界的一个小概览，因此，可能会错过一些你喜爱的美食名人。
我是 数据和植物科学家 在 垂直未来 ，先前在WBC SD。理学硕士。在城市分析和智能城市方面，来自 UCL 的 Bartlett 伦敦的 CASA 和 BSc。在食品科学中从 苏黎世联邦理工学院 。对城市、美食、健康充满热情通过Twitter或LinkedIn联系。
扫一扫
被折叠的 条评论
为什么被折叠?
到【灌水乐园】发言
