通过社交媒体定位自然灾害
每天每时每刻我们都在利用数据。例如，现在当我写这篇文章时，我让网飞流从我观看的节目中检索和收集数据。我看着我在当地克罗格商店买的那瓶金峰减肥茶。为了获得我的燃料积分，我刷了一张卡，记录了我昨天的购物记录，包括那瓶茶。我的苹果 iphone 就坐在我旁边，等着我拿起来记录我的屏幕时间。我的特灵暖通空调恒温器就在我坐的椅子上方，监控着我家的舒适度。即使在我睡觉的时候，我的 fitbit 和 sleep number bed 也会监控我的睡眠模式。
大多数被挖掘的数据是为了个人健康和舒适，或者是企业向我们定制营销其产品的一种方式。数据还可以用于监测和预测地震、森林火灾、恶劣天气等自然灾害。2019 年 3 月 5 日，田纳西州 Maynardville 外仅 4 英里处发生了 3.4 级地震，这是一个我们通常不太听说地震的地区。这次地震发生在被称为南阿巴拉契亚地震带的断层线上，该地震带从阿拉巴马州东北部延伸到田纳西州东部和弗吉尼亚州西南部。这种情况启发我进入 Twitter API，把过去 9 天中所有与包含#地震的推文相关的推文拉下来。在这个故事中，我将提供使用 rStudio 和 R 编程语言下载 Twitter 数据的步骤。
如果你还没有下载 Twitter 数据的必要密钥，第一步就是创建一个账户。为此，请遵循以下步骤。
登录https://developer.twitter.com网站。注册并创建一个帐户。拥有帐户后，请创建一个应用程序。
填写所有细节，创建一个应用程序名称，我用我的谷歌账户【https://myaccount.google.com/的网址。一旦你准确地填写这个，你就会收到你的钥匙。您需要提供一个电子邮件地址，以便 Twitter 可以验证您的地址。
请注意，这是我接收关键信息的地方，我用这些信息从 API 中提取数据。一旦我们得到了我们需要的密钥，我们就可以开始用 r 写代码了。
## install rtweet from CRAN
## install dev version of rtweet from github
devtools::install_github("mkearney/rtweet")# Install the following packages
install.packages(c("ROAuth", "plyr, stringr", "ggplot2","wordcloud"),dependencies=T)# hit enter and this is what you will see.
Restarting R session...> install.packages(c("ROAuth", "plyr, stringr", "ggplot2", "wordcloud"), dependencies = T)
Error in install.packages : Updating loaded packages
> install.packages(c("ROAuth", "plyr, stringr", "ggplot2", "wordcloud"), dependencies = T)
Warning in install.packages :
package ‘plyr, stringr’ is not available (for R version 3.5.1)
trying URL '[https://cran.rstudio.com/bin/macosx/el-capitan/contrib/3.5/ROAuth_0.9.6.tgz'](https://cran.rstudio.com/bin/macosx/el-capitan/contrib/3.5/ROAuth_0.9.6.tgz')
Content type 'application/x-gzip' length 76740 bytes (74 KB)
==================================================
downloaded 74 KBtrying URL '[https://cran.rstudio.com/bin/macosx/el-capitan/contrib/3.5/ggplot2_3.1.0.tgz'](https://cran.rstudio.com/bin/macosx/el-capitan/contrib/3.5/ggplot2_3.1.0.tgz')
Content type 'application/x-gzip' length 3622365 bytes (3.5 MB)
==================================================
downloaded 3.5 MBtrying URL '[https://cran.rstudio.com/bin/macosx/el-capitan/contrib/3.5/wordcloud_2.6.tgz'](https://cran.rstudio.com/bin/macosx/el-capitan/contrib/3.5/wordcloud_2.6.tgz')
Content type 'application/x-gzip' length 229178 bytes (223 KB)
==================================================
downloaded 223 KB
现在我们将加载我们的库。注:当我加载我的库时，我把光标放在我想要加载的每个库之后，然后按回车键。
## load the package libraries
library(rtweet)
library(ROAuth)
library(plyr)
library(ggplot2)
library(wordcloud)
library(tm)# hit enter and you should see this
The downloaded binary packages are in
/var/folders/fc/s97j11jn495g6l0k1ds7x2540000gn/T//RtmpbxuZuS/downloaded_packages
> ## load the package libraries
> library(rtweet)
> library(ROAuth)
> library(plyr)
> library(stringr)
Warning message:
package ‘stringr’ was built under R version 3.5.2
> library(ggplot2)
> library(wordcloud)
Loading required package: RColorBrewer
> library(tm)
Loading required package: NLPAttaching package: ‘NLP’The following object is masked from ‘package:ggplot2’:annotate
接下来我将安装 devtools 包，如果它还没有安装的话。
## install devtools package if it's not already
if (!requireNamespace("devtools", quietly = TRUE))
install.packages("devtools")# You should get something like this below
installing the source package ‘devtools’trying URL '[https://cran.rstudio.com/src/contrib/devtools_2.0.1.tar.gz'](https://cran.rstudio.com/src/contrib/devtools_2.0.1.tar.gz')
Content type 'application/x-gzip' length 388953 bytes (379 KB)
==================================================
downloaded 379 KB* installing *source* package ‘devtools’ ...
** package ‘devtools’ successfully unpacked and MD5 sums checked
** R
** inst
** byte-compile and prepare package for lazy loading
** help
*** installing help indices
*** copying figures
** building package indices
** installing vignettes
** testing if installed package can be loaded
* DONE (devtools)The downloaded source packages are in
‘/private/var/folders/fc/s97j11jn495g6l0k1ds7x2540000gn/T/RtmpbxuZuS/downloaded_packages’
好了，是时候检索 Twitter API 密钥了。我们将把它们放在 R 代码中，以便从 API 中访问和检索数据。
## access token method: create token and save it as an environment
# Note put in your API, not the info below, it will not work.create_token(
app = "chris_kehl",
consumer_key = "dvkGkS6njqwiethinrg74bdt",
consumer_secret = "av6ZVYgrJk8dHPqTrnTwirng",
access_token = "12538373736-sEAAMIHF069nldOEAs3CY87eCEpvt8a4RCD9m",
access_secret = "J5TEGT6rDtenrye9d2O6ZifSjwiAECRp8o9R5x")# hit enter and you will get this> ## access token method: create token and save it as an environment variable
> create_token(
+ app = "chris_kehl",
+ consumer_key = "dvkGkS6njqwiethinrg74bdt",
+ consumer_secret = "av6ZVYgrJk8dHPqTrnTwirng",
+ access_token = "12538373736-sEAAMIHF069nldOEAs3CY87eCEpvt8a4RCD9m",
+ access_secret = "J5TEGT6rDtenrye9d2O6ZifSjwiAECRp8o9R5x")
<Token>
<oauth_endpoint>
request: [https://api.twitter.com/oauth/request_token](https://api.twitter.com/oauth/request_token)
authorize: [https://api.twitter.com/oauth/authenticate](https://api.twitter.com/oauth/authenticate)
access: [https://api.twitter.com/oauth/access_token](https://api.twitter.com/oauth/access_token)
<oauth_app> chris_kehl
key: dvkGkS6njqwiethinrg74bdt
secret: <hidden>
<credentials> oauth_token, oauth_token_secret
是时候检索包含#地震信息的 twitter 数据了。
## search for 18000 tweets using the rstats hashtag
eq.list <- search_tweets(
"#earthquake", n = 18000, include_rts = FALSE
)# output will be something like this> ## search for 18000 tweets using the rstats hashtag
> eq.list <- search_tweets(
+ "#earthquake", n = 18000, include_rts = FALSE
+ )
Downloading [===================================>-----] 88%
>
为了进行分析，我们将对刚刚检索到的数据进行格式化。
# Create a Corpus
eqCorpus <- Corpus(VectorSource(eq.list$text))# set up stemming
eqCorpus <- tm_map(eqCorpus, stemDocument)# output will be something like this> # Create a Corpus
> eqCorpus <- Corpus(VectorSource(eq.list$text))
> # set up stemming
> eqCorpus <- tm_map(eqCorpus, stemDocument)
Warning message:
In tm_map.SimpleCorpus(eqCorpus, stemDocument) :
transformation drops documents
现在让我们来绘制过去九天中收到的推文。
## plot time series of tweets
ts_plot(eq.list, "3 hours") +
ggplot2::theme_minimal() +
ggplot2::theme(plot.title = ggplot2::element_text(face = "bold")) +
ggplot2::labs(
x = NULL, y = NULL,
title = "Frequency of #earthquakes Twitter statuses from past 9 days",
subtitle = "Twitter status (tweet) counts aggregated using three-hour intervals",
caption = "\nSource: Data collected from Twitter's REST API via rtweet"
)# Output will look like this
> ## plot time series of tweets
> ts_plot(eq.list, "3 hours") +
+ ggplot2::theme_minimal() +
+ ggplot2::theme(plot.title = ggplot2::element_text(face = "bold")) +
+ ggplot2::labs(
+ x = NULL, y = NULL,
+ title = "Frequency of #earthquakes Twitter statuses from past 9 days",
+ subtitle = "Twitter status (tweet) counts aggregated using three-hour intervals",
+ caption = "\nSource: Data collected from Twitter's REST API via rtweet"
+ )
我们要分析的情节是:
我们可以分析我们的频率图，看看我们有什么类型的活动。正如所指出的那样，田纳西州在 2019 年 3 月 5 日发生了地震，我们的峰值有所增加，但没有 3 月 1 日那么多。
我们可以搜索检索到的数据，看到一些震级为 5.5、5.3 级的地震活动。通过滚动数据，我们可以看到从阿根廷一直延伸到阿拉斯加的地震。我们现在将绘制最活跃推文的经度和纬度。
## create lat/lng variables using all available tweet and profile geo-location data
eq <- lat_lng(eq.list)## plot state boundaries
par(mar = c(0, 0, 0, 0))
maps::map("world", lwd = .90)## plot lat and lng points onto world map
with(eq, points(lng, lat, pch = 20, cex = .75, col = rgb(0, .3, .7, .75)))
我们运行我们的代码，看到我们的世界图。
从我们的图中，我们可以看到来自加利福尼亚和田纳西的活动。这是大多数 twitter 信息的来源。我们可以查看检索到的数据来分析地震的严重程度。
概括起来
我们可以使用数据来分析任何事情，从谁在网飞上传什么内容，到我们的 iphone 用户每天看屏幕的时间。但是，如果我们可以分析推文来拯救生命，或者向有需要的人提供必要的资源，那又会怎么样呢？在这个故事中，我们通过绘制推特的频率和世界地图上的点来分析地震活动。我们可以通过图中所示的推文频率来分析进一步调查的必要性。我们可以使用我们的示例代码来分析其他灾难，例如叙利亚难民危机、与俄罗斯入侵克里米亚有关的公众情绪，以及监视恶劣天气期间发生的情况。可能性是无穷的，只要把#地震改成#救命！，#克里米亚，#下一场灾难。也许我们可以用这些数据来拯救生命，或者创造改变。请跟随我来到 https://www.medium.com/@chris.kehl 的第二部分。
基于位置的推荐
使用机器学习创建基于位置的推荐器。
Location-based recommendation — Source Microsoft.
推荐系统在不同的应用中被广泛用于预测用户对产品或服务的偏好或评级。在过去的几分钟或几个小时里，你很可能在网上遇到过某种类型的推荐系统或与之互动过。甚至这篇文章可能已经被媒体或其他渠道通过推荐系统推荐给你了。
这些推荐系统可以有不同的类型，其中最突出的包括基于内容的过滤和协同过滤。在本文中，我们将研究基于位置的推荐，其中我们特别关注地理位置，以便利用用户的位置提供更相关的推荐。
为了说明基于位置的推荐器的关键方面，我们将使用 K-Means 算法和来自 Kaggle 的 Yelp 数据集来执行一个简单的基于位置的推荐。这些数据来自 JSON 文件，可以很容易地用 pandas 读取。
下表显示了数据集的前 5 行。在该表中，业务坐标可用于每个业务的评级星和评论数。
Yelp Dataset
探索性数据分析(EDA)和预处理
在本节中，我们将探索和预处理数据集。该数据集包含 yelp 用户的评论，并包括许多类别。为了简化我们对推荐模型的分析和解释，我们将把重点放在餐馆上，您可以随意选择您感兴趣的任何其他类别。
我们首先通过过滤包含单词“餐馆”的所有企业来创建餐馆数据框架。然后，我们绘制餐厅的星级分布图，从 1 到 5。该图清楚地显示，大多数餐厅都有 4.0 或 3.5 星，而只有一小部分餐厅有 5.0 星。
Stars — Countplot
让我们还基于review_count和stars的组合创建一个排序的数据框架，并可视化数据集中的前 20 家餐馆。
这里的选择是不仅根据**stars**对餐馆进行排序，还根据餐馆有多少评论进行排序。这是数据集中排名前 20 的餐厅的图表。
Top 20 restaurants — Fix the plot
最后，在我们进入 K-Means 聚类和基于位置的推荐之前，让我们也探索一下这些餐馆的位置，并想象成一张地图。我们将使用 Plotly Express 库，它提供了一种灵活而简单的方法来绘制带有漂亮底图的地理数据。
正如您在下面的地图中所看到的，该数据集并不完整，仅包含美国和加拿大的几个城市。
Map of the Restaurants in the dataset
为了定制我们的推荐系统，我们将只关注内华达州的拉斯维加斯，它在这个数据集中拥有最多的餐馆。
我们为拉斯维加斯的所有餐馆创建了一个缩放地图。下面你可以看到拉斯维加斯所有餐馆的地图。地图的颜色基于餐厅的星级数，而大小则表示该特定餐厅的评论数。正如你所看到的，高等级的餐馆聚集在市中心。
k 均值聚类
在执行 K-Means 聚类算法之前，让我们首先确定数据集中的最佳聚类数。这里我们使用一种叫做肘法的技术。
肘方法有助于我们确定合适的聚类数。这并不能给出准确的聚类数，但可能有助于选择聚类数。从下图可以看出，它有多个弯头。
Elbow method K-Means Clustering
借助于肘方法，我们选择 5 个聚类，并在坐标上执行 K-Means 聚类算法。我们还在lasVegas数据帧中填充了一个新列。
最后，我们根据 K-Means 算法创建了一个聚类散点图。该图清楚地显示了具有 5 个聚类的算法的输出。
K-Means clustered points
基于位置的推荐
在这一节中，我们将把 K-means 聚类和餐馆评级结合在一起，根据用户的位置推荐餐馆。让我们首先创建一个基于评论数和星级数的 LasVegas 餐馆分类数据框架。我们称之为数据帧top_restaurants_lasVegas.
我们准备创建一个基于用户位置推荐餐馆的功能。我们创建了一个基于用户位置推荐餐馆的功能。
这里，我们的函数将获取纬度和经度以及数据帧。该函数首先预测坐标属于哪个聚类，并将数据帧进一步限制到该聚类，并返回该聚类中的前 5 家餐馆。
这个函数很简单，只需要坐标，但我们可以进一步发展，接受其他输入，如地址或邮政编码。这需要一些地理编码应用程序，所以现在这足以说明基于位置的推荐的威力。在使用我们的推荐函数之前，让我们创建一个测试数据框架。
现在，我们可以使用我们的推荐功能，将数据框架中的任何用户作为输入，来测试我们简单的基于位置的推荐的功能。记住这个测试数据只是为了方便，推荐函数直接接受任何给定的坐标。
对于第一个用户，K-Means 算法预测该聚类，并提取该聚类的前 5 家餐馆并推荐它们。红点是用户位置，蓝点是推荐餐厅。
User 1 Red and Recommended Restaurants Blue
让我们再举一个例子，在这个例子中是用户 3。同样，我们的“推荐餐馆”功能预测该用户所属的聚类，并找出该聚类中的前 5 家餐馆。
结论
在本教程中，我们使用 Yelp 数据集演示了 python 基于位置的推荐。我们用餐馆展示了一个基于位置的推荐应用程序。我们的推荐系统使用基于餐馆坐标的 K-Means 聚类，并且只考虑餐馆的评级(评级数量和星级)来推荐餐馆。我们也可以考虑餐馆的类型来定制推荐。这也可以进一步发展，以纳入其他输入方法，如地址或邮政编码。最后，我们的推荐没有明确考虑用户和推荐餐馆之间的距离，尽管 K-means 聚类限制了推荐餐馆的范围。
本教程的笔记本可以在这个 Github 资源库中找到。
[## shaka som/基于位置的建议
推荐系统在不同的应用中被广泛用于预测用户的偏好或评级…
github.com](https://github.com/shakasom/location-based-recommendations/blob/master/README.md)
位置位置位置
如何使用机器学习和一点黑魔法来创建地理区域嵌入？
Zoopla.co.uk 之旅主要集中在地理区域:这是 Zoopla 搜索的唯一强制输入。消费者在他们的 Zoopla 搜索中引用了许多独特的地理区域，准确地说是去年的 32，237 个**。**
强烈希望更好地了解所有这些领域，并围绕区域发现和区域比较等领域打造更好的产品。但是你如何描述和比较 32237 个地理区域呢？一个地区的哪些特征是相关的，你如何以可靠和可比较的方式收集所有这些数据？不可能的！或者是？
数据科学可以将机器学习和一点黑魔法应用于许多问题。我们的解决方案:让我们直接从 Zoopla 的搜索数据中提取一个区域的相关特征及其值。
我们所需要学习的任何地理区域的描述是一个惊人的简单的数据点。消费者在 Zoopla 搜索中引用的地理区域序列。
zoopla.co.uk property search
使用 Zoopla 搜索的消费者有不同的需求和优先考虑。他们还对不同的地理区域以及他们是否适合租房或买房有所了解。消费者以他们在搜索中参考地理区域的顺序来表达这种知识。搜索序列不是随机的，因此我们可以在机器学习的帮助下提取区域描述。
Zoopla.co.uk area search sequences (pixabay )
我们部署的模型被称为 Skip-Gram ，这是一个试图预测搜索序列的神经网络。该模型使用机器学习来学习区域的描述(称为嵌入)和搜索序列(上下文矩阵)内区域的上下文。
区域描述(嵌入)包含区域的可负担性或交通联系等信息。一个区域的上下文描述了一个区域与其他区域相比如何基于消费者在他们的搜索序列中表达的需求和目标。
综上所述，Skip-Gram 模型根据目标区域描述和消费者搜索环境预测 Zoopla 搜索中的区域顺序，例如年轻的专业人士或有孩子的家庭。
Area Skip-Gram Model
我们用超过 1 亿多个独特的搜索序列来训练神经网络。如下图所示，区域嵌入本身就很有用。此外，嵌入是一种非常有用的方法，可以对分类变量进行编码，以便在其他模型中使用，例如 Zoopla 的推荐引擎。在另一篇博客文章中有更多关于如何构建推荐引擎的内容。
结果:
该模型允许我们根据相似性来比较区域，例如，这是伦敦 out 代码与 SE11 的相似性，用于租赁一处房产:
area similarity to London SE11 for renting a property
这张地图显示了购买房产的类似地区:
Great Britain area similarity to Brighton for buying a property
伯恩茅斯地区与布莱顿非常相似，这并不奇怪:例如，这两个城市都有大型大学。滨海韦斯顿也一样。这三个地区都有通往伦敦的像样的铁路。布拉德福德地区可能是一个更令人惊讶的结果，但布拉德福德也有一所大学，并因靠近利兹和曼彻斯特而获得了额外的相似性。
区域数学的乐趣！
因为机器学习模型给了我们客观的区域数学描述，没有什么能阻止我们做一些区域数学！
搜索哈默史密斯+安吉尔-肯辛顿，我们会看到肖尔迪奇、哈克尼、霍克斯顿和伦敦机场。
另一个例子:有哪些地区与布里克斯顿和卡姆登镇相似，但与克拉彭和 SE1 不同？这种比较是对中产阶级化、支付能力和交通联系的玩弄。
肯特郡是一个直观的结果。西部的费尔瑟姆和特威肯汉姆是正在经历大规模中产阶级化的廉价地区。东边的 Rainham 被认为是梅德韦镇中最负盛名的，像布里克斯顿和卡姆登镇一样有当地的感觉，并且有极好的交通联系。罗姆福德也有很好的联系，并经历了大量的重新开发，预计跨铁路的到来。
信用和来源:
伦敦外部代码形状文件可从位于 https://data.london.gov.uk/的伦敦数据存储中获得，铁路和河流的形状文件可从 OpenStreetMap.org 获得，并可从http://download . geofabrik . de/Europe/great-Britain/England/greater-London . html下载
英国邮政编码形状文件可以从https://www.opendoorlogistics.com/data/和
- 包含皇家邮政数据皇家邮政版权和数据库权利 2015
- 包含国家统计数据皇冠版权和数据库权利 2015。
所有图像都是用 GeoPandas 创建的，并由-NC-SA 4.0 授权 CC，代码如下:
import geopandas as gpd
import matplotlib.pylab as plt
import matplotlib.patheffects as PathEffectssimilarity_scores = compute_cosine_similarities(area_embedding_vectors, se11_embedding_vector)geo_df_outcodes = gdp.read_file("path/to/out-code-shapefile.shp")
geo_df_railways = gdp.read_file("path/to/osm-shape-files/gis_osm_railways_free_1.shp")
geo_df_water = gdp.read_file("path/to/osm-shape-files/gis_osm_water_a_free_1.shp")geo_df_outcodes.join(similarity_scores).plot(
column="score", cmap='Reds', linewidth=0.8, ax=ax, alpha=0.8,
edgecolor='0.5', vmin=0, vmax=1, ax=ax
)
geo_df_railways.plot(alpha=1, edgecolor='0.2', ax=ax, color='black', linewidth=0.1)
geo_df_water.plot(alpha=1, edgecolor='0.5', ax=ax, color='deepskyblue', linewidth=0.8)labels = similarity_scores.nlargest(10, "score").reset_index()['name'].to_list()similarity_scores.apply(
lambda x: ax.annotate(
s=x.name, xy=x.geometry.centroid.coords[0], ha='center'
).set_path_effects([PathEffects.withStroke(linewidth=2,
foreground='w')]
) if x.name in labels else None,
axis=1
)plt.axis('off')
Jan 是公司数据转型方面的成功思想领袖和顾问，拥有将数据科学大规模应用于商业生产的记录。他最近被 dataIQ 评为英国 100 位最具影响力的数据和分析从业者之一。
在 LinkedIn 上连接:【https://www.linkedin.com/in/janteichmann/】
阅读其他文章:https://medium.com/@jan.teichmann
位置，位置，位置！在数据科学中
我们如何接受数据科学中的位置组件？它只是数据集中的另一列吗？或者空间是特殊的。
位置，位置，位置！
你已经听过很多次了。这是房地产界的一个普遍说法。这也适用于数据科学吗？我们如何接受数据科学中的位置组件？它只是数据集中的另一列吗？或者空间是特殊的。
位置数据(大数据)无处不在，因为我们每天都在创建海量的地理数据，从推文和地理标记图像到智能手机跟踪服务。然而，在数据科学中，位置组件并没有深入融入主流。在本文中，我强调了这些松散集成的位置数据科学的一些方面，以及建立紧密集成系统的潜力。
时间序列与位置数据
要意识到我们将位置视为数据科学中的一个补充栏目，只需看看 Kaggle 平台。
为了说明我的观点，让我们来比较两个密切相关的概念:空间和时间。我猜你已经熟悉时间序列分析了。它通过分析一段时间内的序列或预测未来值来提取信息。
即使在空无一物的空间里，时间和空间依然存在。肖恩·卡罗尔
我们都更喜欢数据和事实，所以我在 Kaggle 中搜索了这三个词:“时间序列”、“地理”和“地理空间”。结果一点也不意外。
截至 2019 年 10 月 27 日，在 Kaggle 中搜索时间序列将返回 430 个数据集和 1149 个笔记本。
Time series Kaggle
另一方面，地理搜索产生了几乎相同的数据集编号 391，但看看笔记本和竞赛结果。与时间序列中的 1149 个笔记本相比，包含地理单词的笔记本只有 45 个，如下所示。
你可能会说这些结果不具有普遍性或者有点狭隘，毕竟 Kaggle 只是数据科学项目的一个平台。因此，我查看了谷歌趋势，结果显示了类似的模式——地理数据分析覆盖范围更少，时间序列分析更多，如下所示。
Google Trends
越来越清楚的是，位置数据只是数据集中的另一列，我们通常谨慎地处理它，或者从我们的分析中完全删除它。让我们看看这种超然的原因。
断开位置数据科学
空间数据科学不仅产生漂亮的数据可视化，即地图，而且产生大量的空间思维视角。通过空间分析获得的洞察力为许多应用提供了丰富且无与伦比的洞察力。然而，一些挑战阻碍了将位置组件完全纳入数据科学。
尽管位置数据很重要，但我们今天看到数据科学和地理空间组件之间的融合较少。这是为什么呢？
除了一些例外，空间部分在主流数据科学中大多被忽略。一个很好的例子是，位置组件是如何嵌入深度学习应用程序中的，而无需关注空间方面。即使是卫星图像中的大多数深度学习应用，一旦我们将它们输入神经网络，也会失去空间部分。
此外，位置分析和数据科学有很大的重叠，但在很长一段时间内没有联系。一个早期的例子是约翰·斯诺的霍乱分析。随着 1854 年伦敦霍乱的爆发，约翰·斯诺绘制了死亡的地理位置图，并在两个水源发现了一群死亡病例。可以说，这是最早记录的数据科学和数据分析案例之一。
空间数据科学教授吕克·安瑟林最近将这种脱节称为**“空间怀疑主义”。**他反映说，大数据和机器学习社区中的许多人仍然不相信将空间视角紧密集成到数据科学中的重要性，就像 John snow 的霍乱分析一样。
另一个主要障碍是由于历史文物。处理地理数据在计算上很困难，但现在不再是这样了。如今，空间数据库的主流集成以及操作大地理信息的工具已经出现。例如，谷歌的 BigQuery 就包含了 GIS 功能。
此外，很少有应用程序和项目从简单的位置整合到数据科学中，转向完全集成的地理数据科学。Spatial Hadoop 和最近发布的Cu spatial——使用 NVIDIA RAPIDS DataFrame 库在 GPU 上加速——只是一些进步。
走向地理数据科学
地理数据科学是明确关注数据科学的位置(空间)部分的学科。它带来了数据科学领域中特定于地理数据的理论、概念和应用。我们都知道东西在哪里以及如何进入它们的重要性，这影响到从商业物流、运输、气候变化到与选址决策相关的一切。
我认为是时候超越数据科学和位置组件之间的脱节了。我们需要认识到，将数据可视化为地图并不等同于空间数据科学。正如贾维尔·德·拉·托雷上周在空间数据科学会议上所说，是时候从在地图上可视化数据转向使用地图分析数据了。
我对空间数据科学充满热情，如果你喜欢讨论空间数据科学，可以在 Twitter 上通过@shakasom 与我联系。
黄土
使用局部回归平滑数据
Photo by Vinícius Henrique on Unsplash
如果你对一个物理现象产生的数据进行采样，你会得到噪音。测量信号的传感器会将噪声添加到信号中，或者噪声是产生数据的过程固有的随机性。最近，我不得不处理一个由汽车引擎产生的噪音数据流，并需要找出一种方法来过滤噪音。由于信号产生过程的物理性质，采样频率不是恒定的，因此排除了任何基于频率的噪声滤波技术。我需要找到一种方法来滤除噪声，并重建信号以供进一步处理。
为了从测量的噪声中恢复信号，我们必须首先对噪声是如何产生的做一些假设。从统计学的角度来说，这意味着我们必须假设噪声的某种分布，即噪声产生的数学描述。最常见的假设包括根据高斯分布产生的随机噪声、将噪声添加到信号中的加性模型，以及独立于 x 的误差项，如下所示:
Additive noise
顺便提一下，这与线性回归模型通常假设的噪声模型完全相同。这里的噪声模型是这样的:
The linear regression noise model
从某种意义上说，通过将模型拟合到您的数据，您正在尝试从通过数据方差实现的测量噪声中恢复潜在信号。
那么，我应该用线性回归来平滑信号吗？如果信号本身可以用一个线性函数来建模，这个函数可能包含非线性项，那么答案是肯定的。但在这个特定的场景中，我将处理一个高度非线性的信号，它反映了一个配送工具的日常操作:大量的平坦周期夹杂着可变长度的纯粹疯狂周期。用单一的线性回归模型处理整个信号是不可能的。但是…
如果我不使用同一模型处理整个信号，而是使用同一类型的不同模型来平滑信号的小区域和局部区域，会怎么样？为什么不把这个想法进一步发展，为我们需要平滑的每个点考虑一个专门的模型呢？
考虑当地情况
看我的锦囊妙计，发现一个老朋友:黄土 —局部加权运行线平滑器。这是非参数平滑器，尽管它的核心使用线性回归。与任何平滑器一样，该算法的思想是从噪声样本中恢复固有信号。
那么黄土是怎么起作用的呢？让我们从一个像下面这样的噪声信号开始。
Noisy signal
这是一个合成生成的正弦波，带有附加的高斯噪声。正弦波以红色绘制，而噪声样本显示为蓝点。为了模拟不规则采样的信号，从均匀分布中随机采样 x 值，并适当缩放。使用添加了高斯噪声的正弦函数来计算相应的 y 值。
那么我们如何从蓝点到红线的近似值呢？首先，将红线视为一个有序的等间距 x 值序列，在本例中位于 0 和 2π 之间。对于这些值中的每一个，选择采样点的适当邻域，并将它们用作线性回归问题的训练集。使用生成的模型，估计您的点的新值。让我们更详细地探讨一下这个想法
黄土背后的理念
该算法以逐点方式估计潜在函数。对于 x 的每一个值，我们通过使用其相邻的采样(已知)值来估计 f(x) 的值。这非常类似于 KNN 算法，其中 k 、窗口大小是可调参数，并且在这种特定情况下，将确定所得估计的平滑度。从某种意义上来说， k 就是你的偏向 vs 。方差旋钮。 k 的大值将导致更高的偏差，而低值将导致更高的方差。
第一步是收集我们想要估算的的值 x 。姑且称这些x’和y’。通过给黄土算法输入**【x’，并使用采样的 x 和y 值，我们将获得一个估计值*【y’*。从这个意义上说，黄土是一种非参数算法，必须使用所有数据集进行估计。
既然我们有了x’，我们必须使用简单的欧几里德距离找到它的 k 最近邻。让我们把产生的有序集称为 D 。
下一步将集合 k 距离的集合 D 转换为有序集合 W ，该集合包含稍后将在线性回归过程中使用的权重。使用专门的权重函数来计算这些权重，该权重函数根据其到x’的距离来为的每个 k 邻居分配重要性。
权重函数
使用三次函数计算距离权重:
Tri-cubic weighting function
该函数看起来像一顶帽子，仅在 -1 和 1 之间有正值。在这个区间之外，函数为零。该函数如下所示:
Graph of the tri-cubic weight function
由于该函数仅对-1<x<1 有正的结果，我们必须通过将距离除以在 D 中观察到的最大值来归一化距离。更具体地说，
Weighting function
在这里，我们将 d(x，x’) 表示为xk最近邻之一与 x’ 之间的距离。标准化的效果是距离越大，权重越低。在极端情况下，与最大距离对应的点的权重为零，距离为零的点的权重可能最高，为 1。这就是“ locality ”效果是如何实现的，方法是为最接近我们想要计算预测的位置的训练数据分配更高的重要性。
顺便提一下，你可能会发现这个函数与三次核函数有惊人的相似之处。这些函数之间的规模差异(70/81)与核函数必须在其域上积分的要求有关，而这里的要求是宽松的。
我们现在准备使用简单的加权线性回归来计算估计值，该回归使用来自 D 的 x 值以及相应的 y 值进行训练。
线性回归
线性回归是监督机器学习方法的基础。很有可能，你开始了你的 ML 之旅，学习这种方法的内部原理，考虑到波特兰家庭的物理特征，你可能试图计算出他们的销售价格。或者完全是别的原因，但你知道该怎么做，不是吗？
碰巧的是，线性回归的一个特殊版本，加权线性回归，是黄土的核心。对于我们开始估计的每个点(【x’)，黄土算法必须建立一个线性回归模型，该模型将计算相应的输出(【y’)，使用*【k】******【x’*的最近邻居和一组对其重要性进行评级的权重。
局部线性回归通常模拟低维多项式、直线或二次多项式。
The first-degree regression equation
The second-degree regression equation
加权线性回归是一个已知的问题，在网上有大量的记录。由于要解决的问题的典型低维数，我们将求助于闭合形式的正规方程进行参数估计。在未加权的情况下，这些方程是:
Normal equations for linear regression
Were beta 是线性参数的向量，X 是包含所有 X 个观测值的矩阵，排列如下:
The X matrix
具体地说，这个矩阵用 n 维数和 m 观测值对样本建模。请注意，我在矩阵的第一列中包含了截距项。对于我们对二次多项式建模的情况，这个矩阵实际上是:
Second-degree X matrix
一旦我们有了β向量，就可以使用下面的等式来计算新的值 y :
将这一概念扩展到使用权重实际上非常简单，正常方程只需要一个额外的项:
Weighted normal equation
这里，权重矩阵 W 在对角线上具有所有计算的权重，所有其他元素被设置为零。不幸的是，正如您将在实现的 Python 代码中看到的，矩阵方法可能有点慢。如果您坚持一级模型，可以使用更简单的数学方法采取替代方法:
****
看起来很复杂，但通过使用向量的内积来消除显式求和，实现要简单得多。批注注: d 代表 D 中的项数，实际上是 k 。
继续编码。
Python 库
您可以在 StatsModels Python 包中找到这个平滑器的实现。通读方法文档，您会看到 lowess 函数返回一个与两个输入数组( x 和 y )维数相同的数组。这意味着只有观察值是平滑的，所以如果你需要中间的任何其他值，你将不得不以某种方式对它们进行插值。但是事情不应该是这样的。
Python 实现
对于本文，我开发了一个基于 NumPy 的新实现，它利用了它的矢量化特性，代码可以在这个 GitHub 库中找到。开发代码时考虑了矢量化，函数中只有一个循环来确定最接近值的索引。让我们浏览一下代码，看看它是如何工作的。
三次加权函数是完全矢量化的，它处理 x 值的数组。首先，输出数组 y 被创建为与输入数组 x 具有相同的维数。接下来，创建一个索引数组来加强函数的定义域，最后，计算函数本身。请注意，索引数组在输入和输出数组中都使用。我的第一种方法是使用 Numba 对代码进行矢量化，但是后来我意识到这种方法具有相同的性能，并且消除了不必要的编译。
初始化时，两个输入数组都必须归一化，以避免失去重要性(也称为灾难性消除)的问题。这是通过重新调整到 0 和 1 之间的间隔来完成的。
在包含最小距离窗口的索引的索引数组的帮助下，从距离数组计算权重。
该索引数组在下一个函数中计算:
为了计算到*【x’的最小总距离的范围，我们首先确定距离数组中最小距离的索引。知道索引必须是连续的，我们可以使用这个初始索引作为一个不断增长的索引列表的根。当最小索引位于距离数组的极值时，函数顶部的测试只处理边缘情况。下面的循环增加了索引列表，从最小距离的索引开始，根据需要在左侧和右侧添加项目，并保持列表自然排序，插入到左侧并追加到右侧。注意循环次数限制为kT7-1*。****
现在，我们进入代码的核心。估算 f(x)的函数可用于两种模式:矩阵模式或统计模式。在矩阵模式下，可以指定多项式次数，但性能会较低。统计代码更快，但只能模拟线条。
该功能从标准化输入 x 值开始，并计算其与所有训练值的距离。距离数组与训练数据具有相同的维数。接下来，找到最小距离范围并计算相应的权重。注意权重数组有 k (窗口大小)项。最后，训练回归，并使用上述任一方法计算的估计值。请注意，如果您想要使用多项式回归，代码将使用“矩阵模式”。
最后，这里有一个如何使用代码的示例(数据值取自 NIST ):
请注意，除了训练数据中的值，您还可以提供 x 的值。
以下是对与第一个图表相同的数据使用平滑函数的示例:
Interpolated function in black
您可以使用 GitHub repo 中的配套笔记本来玩这个图表。
结论
我们已经讨论了使用黄土局部回归模型的基本原理，并揭开了它如何工作的面纱。开发并展示了一个 Python 实现，它大量使用了 NumPy 库及其矢量化特性。请自行使用 GitHub 资源库中的代码，并在评论中告诉我您的想法。
感谢您的宝贵时间！
笔记
- 我在[1]中找到了这个定义。作者没有提到“LOWESS”一词。
参考
[1] Gareth，J. Witten，D. Hastie，T. Tibshirani，R. (2013 年)。统计学习介绍及在 R 中的应用。纽约:施普林格
[2]阿尔佩登，E. (2014 年)。机器学习入门。第三版。麻省剑桥:麻省理工学院出版社。
[3]斯塔默，J. (2017)。 StatQuest:对数据拟合曲线，又名 lowess，又名黄土， YouTube 。
** [## joo Paulo Figueira-数据科学家- tb.lx by 戴姆勒卡车和公共汽车| LinkedIn
查看 joo Paulo Figueira 在全球最大的职业社区 LinkedIn 上的个人资料。圣保罗列出了 1 份工作…
www.linkedin.com](https://www.linkedin.com/in/joao-paulo-figueira/)**
日志——K-均值聚类的距离测量方法指南
在本指南中，我试图介绍 K-Means 聚类中可以使用的不同类型和特征的距离
让我们从集群的简单介绍开始。聚类是将数据点分成若干组的任务，使得同一组中的数据点比其他组中的数据点更类似于同一组中的其他数据点。简而言之，目标是隔离具有相似特征的群体，并将他们分配到集群中。
K-Means 聚类是众多聚类算法中的一种。背后的想法是定义集群，使集群内的总变化(称为集群内总变化)最小化。 K-means 算法可以总结如下:
1\. Specify the number of clusters(k) to be created.2\. Select randomly k objects from the data-set as the initial cluster centers or means.3\. Assign each observation to their closest centroid, based on the specified distance[the type of distance is what we will be exploring in this article, in the above case it is Euclidean] between the object and the centroid.4\. For each of the k clusters update the *cluster centroid* by calculating the new mean values of all the data points in the cluster. The centroid of a *K-th* cluster is a vector of length *p* containing the means of all variables for the observations in the *K-th* cluster; *p* is the number of variables.5\. Iteratively minimize the total within sum of square. That is, iterate steps 3 and 4 until the cluster assignments stop changing or the maximum number of iterations is reached.
K — Means Clustering visualization [source]
在 R 中，我们通过下式计算 K 均值聚类:
Kmeans(x, centers, iter.max = 10, nstart = 1, method = "euclidean")where
x > Data frame
centers > Number of clusters
iter.max > The maximum number of iterations allowed
nstart > How many random sets of center should be chosen
method > The distance measure to be usedThere are other options too of calculating kmeans clustering but this is the usual pattern.
有不同的距离计算方法像 欧几里德、最大值(切比雪夫距离)、曼哈顿、汉明、堪培拉、皮尔逊、abspearson、abscorrelation、spearman 或 kendall *。*那么如何选择使用哪一个呢？
数据点之间的距离测量
这些方法分为两组，一组基于捕获几何间距，另一组依赖于**相关性。**我们将逐一了解。
G 几何分离
欧几里得，曼哈顿&最大值(切比雪夫)距离
首先，这一部分的很多材料都是从divingintodatascience 现在的离线页面 中引用的，这个网站已经帮了大忙了。
闵可夫斯基距离是一种度量，它告诉我们空间中两点之间的距离。现在闵可夫斯基距离有不同的顺序，我们将很快看到它的含义，我们也将看到为什么我谈论它，而不是欧几里德和其他距离。
2 点 p 和 q 的闵可夫斯基距离的通用公式:
由下式给出:
Minkowski distance
闵可夫斯基距离通常与 r 一起使用，r 为 1 或 2，分别对应于曼哈顿距离和欧几里德距离。在 r 达到无穷大的极限情况下，我们得到切比雪夫距离。
Euclidean distance
Manhattan distance
Maximum(Chebychev) distance
更容易理解的方法是下图
Euclidean(green) vs Manhattan(red)
曼哈顿距离通过聚合每个变量之间的成对绝对差来捕捉两点之间的距离，而欧几里德距离通过聚合每个变量中的平方差来捕捉两点之间的距离。因此，如果两个点在大多数变量上接近，但在其中一个变量上差异更大，欧几里德距离将夸大这种差异，而曼哈顿距离将对此不屑一顾，因为它受其他变量的接近程度的影响更大。切比雪夫距离计算一对数据点特征之间的最大绝对差异。
曼哈顿距离应该给出更稳健的结果，而欧几里德距离可能会受到异常值的影响。这同样适用于闵可夫斯基距离公式中“p”的较高值。随着 p 值的增加，距离度量变得更容易失去稳健性，少数维度中的异常值开始支配距离值。
如果我们使用这些不同的距离度量而不是默认的欧几里德距离度量来画一个**‘圆’，可以观察到它们之间的差异。我们知道圆**是一个与给定点(圆心)等距的点的轨迹。现在，如果我们使用曼哈顿或切比雪夫距离度量来测量点与中心的距离，我们实际上得到的是“正方形”而不是通常的“圆形”圆。
堪培拉距离
这是曼哈顿距离的加权版本。它测量一对数据点特征之间的绝对分数差之和，当两个坐标最接近零时，它对微小的变化非常敏感。
Canberra distance
海明距离
对于分类变量(男性/女性，或小型/中型/大型)，如果两个点属于同一类别，我们可以将距离定义为 0，否则为 1。如果所有的变量都是分类的，那么您可以使用汉明距离来计算错配的数量。
您还可以将分类变量扩展为指示变量，每个变量对应一个级别。
如果类别是有序的(如小/中/大),使得一些类别比其他类别彼此“更接近”,那么您可以将它们转换成数字序列。例如，(小/中/大)可能映射到(1/2/3)。然后，您可以使用欧几里得距离或其他距离来表示定量数据。
马氏距离
我们可以将从一个点到其相应聚类中心的 Mahalanobis 距离视为其欧几里得距离除以该点方向上的方差的平方根。Mahalanobis 距离度量优于 Euclidean 距离度量，因为它允许聚类的结构具有一定的灵活性，并且考虑了变量之间的方差和协方差。
当您使用欧几里德距离时，您假设分类具有同一性协方差。在 2D，这意味着你的集群是圆形的。显然，如果数据中自然分组的协方差不是单位矩阵，例如在 2D，聚类具有椭圆形状的协方差，那么在欧几里得上使用马氏性将是更好的建模。
Mahalanobis distance for a two dimensional vector with no covariance
基于相关性的距离
如果两个对象的特征高度相关，则基于相关性的距离认为这两个对象是相似的，即使观察到的值在几何距离上可能相距很远。当两个物体完全相关时，它们之间的距离为 0。如果您想要识别具有相同总体轮廓的观察集群，而不管它们的大小，那么您应该使用基于相关性的距离作为相异度量。
如果选择了欧氏距离，则具有高特征值的观测值将被聚类在一起。这同样适用于具有低特征值的观测值。
皮尔逊相关距离
Pearson correlation 测量两个简档之间的线性关系的程度。皮尔逊相关分析法是最常用的方法。它也称为参数相关性，取决于数据的分布。该距离基于皮尔逊相关系数，该系数是根据样本值及其标准差计算的。相关系数 ’ r '取值范围从–1(大的负相关)到+1(大的正相关)。
Pearson correlation distance
这一距离还有其他一些变化:
- **绝对皮尔逊相关距离:**在这个距离中，使用皮尔逊相关系数的绝对值；因此，相应的距离位于 0 和 1 之间。
- 非中心相关距离:除了在非中心相关的表达式中样本均值被设置为零之外，这与皮尔逊相关相同。不居中相关系数在–1 和+1 之间；因此距离在 0 和 2 之间。
- **绝对无中心相关距离:**这与绝对皮尔逊相关相同，除了在无中心相关的表达式中样本均值被设置为零。不居中相关系数在 0 和+1 之间；因此，距离在 0 和 1 之间。
艾森余弦相关距离
这是皮尔逊相关性的一个特例，其中 x 和 y 都被零取代:
斯皮尔曼&肯德尔相关距离
**两个变量之间的 Spearman 相关性等于这两个变量的等级值之间的 Pearson 相关性；皮尔逊相关性评估线性关系，而斯皮尔曼相关性评估单调关系(无论是否线性)。**如果没有重复的数据值，当每个变量都是另一个变量的完美单调函数时，就会出现+1 或-1 的完美 Spearman 相关性。
直观地说，当观察值具有相似的(或相同的相关性为 1) 等级(即观察值在变量内的相对位置标签:第一、第二、第三等)时，两个变量之间的 Spearman 相关性会很高。)在两个变量之间，当观察值在两个变量之间具有不同的(或者对于 1 的相关性完全相反的)等级时为低。
Kendall tau 等级距离是计算两个等级列表之间成对不一致的数量的度量。距离越大，两个列表越不相似。Kendall tau distance 也称为冒泡排序距离，因为它等于冒泡排序算法将一个列表与另一个列表按相同顺序放置所需要的交换次数。
斯皮尔曼系数对连续和离散序数变量都适用。斯皮尔曼的ρ和肯德尔的τ都可以表述为一个更 一般相关系数 的特例。
关于 k-均值聚类的几点建议
- 距离测量的值与测量的尺度密切相关。因此，变量通常在测量观察间差异之前进行缩放。当变量以不同的标度测量时(例如:千克、千米、厘米……)，这是特别推荐的；否则，所获得的相异度将受到严重影响。**标准化使得四种距离测量方法——欧几里得、曼哈顿、相关和艾森——比使用非转换数据更相似。**注意，当数据被标准化后，皮尔逊相关系数&欧氏距离之间存在函数关系。
- **k-means 适用于连续变量。不应该使用混合类型的数据。当您的数据包含混合类型的变量时，您可以尝试使用高尔距离。**这里有高尔距离的概述。
没有最佳距离度量。
对于给定的数据集，只有最佳距离度量。距离度量的选择将影响您的聚类，但这取决于数据集和目标，哪种距离度量最适合您的特定应用。
Similarity Measures for continuous data
参考
- https://pdfs . semantic scholar . org/B3 B4/445 CB 9a 2 a 55 fa 5d 30 a 47099335 B3 f 4d 85 DFB . pdf
- https://www . data novia . com/en/lessons/clustering-distance-measures/
- https://stats . stack exchange . com/questions/81481/why-does-k-means-clustering-algorithm-use-only-euclidean-distance-metric
- https://arxiv.org/ftp/arxiv/papers/1405/1405.7471.pdf
- https://stats . stack exchange . com/questions/130974/how-to-use-both-binary-and-continuous-variables-together-in-clustering
- 维基百科(一个基于 wiki 技术的多语言的百科全书协作计划ˌ也是一部用不同语言写成的网络百科全书ˌ 其目标及宗旨是为全人类提供自由的百科全书)ˌ开放性的百科全书
- http://www.divingintodatascience.com/
日志——通过 Python 实现 Excel 和 Outlook 电子邮件交付自动化指南
本文分为两个部分，第一部分处理从 excel 生成图像/PDF，下一部分将它附加到 outlook 电子邮件中并发送出去。它在宏和处理不同的邮箱方面也有额外的好处。
在我们的日常活动中，我们经常会遇到发送大量邮件的任务，这个项目就是这种需求的产物。
首先，我们将记下这个自动化过程所需的步骤，然后尝试自动化其中的每个步骤。这种自动化所需的步骤如下:
1\. Generate image/pdf from an excel.
2\. Create an email with the image/pdf as an attachment.
3\. Send out the emails from outlook.
这看起来很简单，让我们开始吧。从现在开始，试着把重点放在代码上，我会在最后附上一个 GitHub 要点。
步骤 1 —图像/PDF 创建
我们有一个 excel 包含了 GOT 所有战役的数据，还有一个简单的图表显示了各地区的战役数量。在这一步中，我们的目标是从 excel 中创建下图的 PDF/图像。你可以从这里下载 excel。
Graph which we will copy as image
首先，我们将连接到 Excel 应用程序并加载我们的工作表，Python 有许多选项来本地创建常见的 Microsoft Office 文件类型，包括 Excel、Word 和 PowerPoint。然而，在某些情况下，使用纯 python 方法来解决问题可能太困难了。幸运的是，python 拥有名为 pywin32 的“Python for Windows Extensions”包，可以让我们轻松访问 Window 的组件对象模型(COM)并通过 Python 控制微软应用。
xlApp = win32.DispatchEx('Excel.Application')
wb = xlApp.Workbooks.Open(file_location + filename)
ws = wb.Worksheets('Graph') **# name of the sheet in which graph is present.**
完成后，我们将给出想要转换为图像的工作表部分的单元格引用。图像功能由 PIL 包提供。以下代码将单元格引用作为图像复制到剪贴板并保存。您可以使用“质量”参数调整图像质量，请注意，图像质量仅适用于 JPEG 类型的文件。默认图像质量设置为 75，当您增加或减少该值时，图像的大小也会改变。具体可以参考 PIL 图书馆的官方文档。
win32c = win32.constants
ws.Range("B2:I16").CopyPicture(Format=win32c.xlBitmap) **# give the cells for which you need the image**
img = ImageGrab.grabclipboard()
img.save(file_location + 'Graph.jpeg',quality=55) **# image quality and file size is directly linked**
到目前为止，我们已经保存了我们的形象。现在，如果你想将图像转换成 pdf 格式，你可以使用 img2pdf 库。下面的代码读取上面保存的图像，并将其保存为 PDF 格式。如果你需要的只是一张图片，你可以跳过这一步。
**#storing pdf path**
pdf_path = file_location + "Graph.pdf"**#opening image**
image = Image.open(image_path)**#converting into chunks using img2pdf**
pdf_bytes = img2pdf.convert(image.filename)**#opening or creating pdf file**
file = open(pdf_path, "wb")**#writing pdf files with chunks**
file.write(pdf_bytes)**#closing image file**
image.close()**#closing pdf file**
file.close()
完成这一步后，我们有一个图像和一个 PDF:
Image & PDF created
我们理想地完成了第一步，但我想补充一点小奖励。通常，我们有一个基于宏的 excel，而不是常规的 excel，需要在拍摄图像之前运行宏来刷新数据，这也是可能的。您只需要知道要运行的宏的名称。下面的代码将打开 excel 运行指定的宏并保存 excel。
xlApp = win32.DispatchEx('Excel.Application')
wb = xlApp.Workbooks.Open(file_location + filename)**# Give the proper macro name**
xlApp.Application.Run("Module.Macroname")
wb.Save()
wb.Close()
xlApp.Application.Quit()
步骤 2—创建电子邮件
我们将使用与 excel 相同的 pywin32 包来连接 Outlook。
outlook = win32.Dispatch('outlook.application')
T4:这是另一笔奖金。当配置了多个电子邮件 id 时，通常需要从特定的 outlook 邮箱发送邮件，我们将在下一步处理这个问题。如果您想从默认邮箱发送邮件，可以跳过这一步。以下代码将遍历您所有已配置的帐户，当它找到“abc@mail.com”时，将从那里发送邮件。但是请注意，如果它没有获得指定的邮件 id，它将只从默认 id 发送邮件。
sendfromAC = *None
for* oacc *in* outlook.Session.Accounts:
**#give the mail id from which you want to send**
*if* oacc.SmtpAddress == "abc@mail.com":
sendfromAC = oacc
*break* ##----------------------------------------------------------------
mail = outlook.CreateItem(0)
*if* sendfromAC:
mail._oleobj_.Invoke(*(64209, 0, 8, 0, sendfromAC))
##----------------------------------------------------------------
我们的 outlook 对象已经创建，现在我们将添加到/CC / 主题和附件。对于附件部分，我将展示 2 件事情，将 PDF 文件作为文件附在邮件中，将图片嵌入邮件正文中。这应该涵盖 2 个用例。
mail.To = 'abc@mail.com'
mail.Cc = 'abc@mail.com;bcd@mail.com'
mail.Subject = 'Demo Outlook Automation mail'
附上 PDF:
mail.Attachments.Add(file_location + "Graph.pdf")
嵌入邮件正文:
为了在邮件正文中嵌入图像，我们首先需要将它作为普通附件附加，然后将其作为 HTML 邮件正文的一部分嵌入。以下部分处理附加图像和设置属性，以便 outlook 可以正确处理图像。我们很快会谈到嵌入部分。
attachment1 = mail.Attachments.Add(file_location + 'Graph.jpeg')
attachment1.PropertyAccessor.SetProperty("http://schemas.microsoft.com/mapi/proptag/0x3712001F", "Graph")
mail.HTMLBody = "<HTML lang='en' xmlns='http://www.w3.org/1999/xhtml' xmlns:o='urn:schemas-microsoft-com:office:office'> " \
+ "<head>" \
+ "<!--[if gte mso 9]><xml> \
<o:OfficeDocumentSettings> \
<o:Allowjpeg/> \
<o:PixelsPerInch>96</o:PixelsPerInch> \
</o:OfficeDocumentSettings> \
</xml> \
<![endif]-->" \
+ "</head>" \
+ "<BODY>"
接下来，我们必须创建电子邮件正文。您可以将邮件正文设置为普通文本格式或 HTML 格式。我个人更喜欢 HTML 正文格式，因为它给了我很多控制电子邮件应该如何显示，而不是让 Outlook 来判断。下面的代码非常明显，它将设置邮件正文，最后的 ‘cid:Graph’ 部分将在邮件正文中嵌入附加的图像。
mail.HTMLBody = mail.HTMLBody + "<BR>Hello,<b> </b>" \
+ "<BR><BR> Refer to the image below: </b> ."\
+ "<html><body><img src='cid:Graph'></body></html>"
步骤 3 —发送邮件
现在我们已经准备好了邮件对象，剩下的就是发送邮件了。这将使用简单的命令来完成:
mail.Send()
万岁！！
Sent Mail
完整的 GitHub 要点补充如下:
已知问题
问题 1
很少会出现脚本无法连接到 excel 的问题。错误是这样的:
module 'win32com.gen_py.00020813-0000-0000-C000-000000000046x0x1x9' has no attribute 'CLSIDToClassMap'
解决方案是从这个堆栈溢出 URL 中派生出来的
这些步骤是:
1\. **Delete the entire gen_py folder**, in my case it was present in this path: C:\Users\'user_name'\AppData\Local\Temp\gen_py folder2\. **Rerun the script**, the above caches will be regeneratedResolving step 2 often results in Issue 2, which has to be resolved in turn.
问题 2
这个问题比问题 1 更常见，运行几天后，代码会抛出一个错误:
File "C:/Users/'user_name'/PycharmProjects/Mailer/mailDispatcher.py", line 64, in MailDispatcher ws.Range("A1:AF25").CopyPicture(Format=win32c.xlBitmap) File "C:\Users\'user_name'\PycharmProjects\mail\venv\lib\site-packages\win32com\client\__init__.py", line 178, in __getattr__ raise **AttributeError(a)**
**in the** *xlApp = win32.DispatchEx('Excel.Application')* **line.**
这个解决方案是通过一种试凑法得出的。这些步骤是:
1\. Replace **xlApp = win32.DispatchEx('Excel.Application')** with **xlApp = win32.gencache.EnsureDispatch('Excel.Application')** *[This is already added in code.]*2\. Run the code for once, this will result **in a mail with corrupted attachments.**3\. **Undo changes of Step 1.**4\. Run the code again, now everything should work.
我真的不知道是什么导致了这些问题，广泛的谷歌搜索没有任何有价值的东西。看起来像是某种缓存损坏。
参考资料:
https://www . ka ggle . com/mylesoneill/game-of-thrones # battles . CSV
https://pbpython.com/windows-com.html
太多堆栈溢出问题
日志——假设检验指南
这是假设检验指南。我试图用一个循序渐进的例子来涵盖理论和实际实现的基础。
对于我们所做的任何研究，我们基本上都是在试图回答一个问题或假设。回答这种问题/假设的一种方法叫做假设检验或显著性检验。
假设检验的结构
假设检验的作用是提供一条路径，我们可以沿着这条路径有效地检验或证明我们的问题/假设。以下是进行假设检验时通常遵循的步骤:
1\. Define the **research hypothesis** for the study.2\. Explain how you are going to **operationalize** (that is, **measure** or **operationally define**) what you are studying and set out the **variables** to be studied.3\. Set out the **null** and **alternative hypothesis** (or more than one hypothesis; in other words, a number of hypotheses).4\. Set the **significance level**.5\. Make a **one** or **two-tailed prediction**.6\. **Select** an appropriate **statistical test** based on the variables you have defined and whether the distribution is normal or not.7\. **Run** the **statistical tests** on your data and **interpret** the **output**.8\. **Reject** or **fail to reject** the **null hypothesis**.[Source](https://statistics.laerd.com/statistical-guides/hypothesis-testing-2.php)
可能会有一些变化，但在大多数情况下，这是遵循的结构。
一个例子
研究假设
一家公司需要为他们的销售人员购买公司汽车。他们联系了当地的福特经销商，询问 2007 款福特探险者 4X4 XLT 的价格。经销商表示，这款车在马里兰州的平均价格(μ)为 29，995 美元。现在，在购买之前，该公司想确认这是不是真的。
使用于操作
该公司对马里兰州的 20 家福特经销商进行了随机抽样，发现平均价格为 30，474.80 美元，样本的标准偏差为 1，972.59 美元。
+------------------------+---+--------+
| [www.academyford.com](http://www.academyford.com) | | 27,595 |
+------------------------+---+--------+
| [www.advantageford.com](http://www.advantageford.com) | | 30,250 |
+------------------------+---+--------+
| [www.appleford.com](http://www.appleford.com) | | 32,705 |
+------------------------+---+--------+
| [www.bobbellford.com](http://www.bobbellford.com) | | 28,485 |
+------------------------+---+--------+
| [www.crouseford.com](http://www.crouseford.com) | | 31,295 |
+------------------------+---+--------+
| [www.darcars.com](http://www.darcars.com) | | 30,075 |
+------------------------+---+--------+
| [www.diehlsford.com](http://www.diehlsford.com) | | 30,505 |
+------------------------+---+--------+
| [www.hagertownford.com](http://www.hagertownford.com) | | 30,820 |
+------------------------+---+--------+
| [www.hillsford.com](http://www.hillsford.com) | | 25,995 |
+------------------------+---+--------+
| [www.hindersford.com](http://www.hindersford.com) | | 32,830 |
+------------------------+---+--------+
| [www.huntford.com](http://www.huntford.com) | | 31,875 |
+------------------------+---+--------+
| [www.koonsford.com](http://www.koonsford.com) | | 32,285 |
+------------------------+---+--------+
| [www.norrisford.com](http://www.norrisford.com) | | 32,580 |
+------------------------+---+--------+
| [www.pittsvilleford.com](http://www.pittsvilleford.com) | | 28,915 |
+------------------------+---+--------+
| [www.plazaford.com](http://www.plazaford.com) | | 29,940 |
+------------------------+---+--------+
| [www.prestonford.com](http://www.prestonford.com) | | 31,720 |
+------------------------+---+--------+
| [www.ramseyford.com](http://www.ramseyford.com) | | 30,265 |
+------------------------+---+--------+
| [www.shafferford.com](http://www.shafferford.com) | | 32,555 |
+------------------------+---+--------+
| [www.towsonford.com](http://www.towsonford.com) | | 31,580 |
+------------------------+---+--------+
| [www.watertownford.com](http://www.watertownford.com) | | 27,226 |
+------------------------+---+--------+
无效假设和替代假设
为了进行假设检验，我们需要将我们的研究假设表述为无效的替代假设。零假设和替代假设是关于人群中发生的差异或影响的陈述。我们将使用我们的样本来测试哪个陈述(即零假设或替代假设)最有可能(尽管从技术上来说，我们针对零假设来测试证据)。
零假设本质上是“魔鬼代言人”的立场。也就是说，它假设我们试图证明的任何事情都没有发生。
零假设(₀ ):
The null hypothesis states that our company believes that the average cost of a 2007 Ford explorer is more than or equal to the quoted $29,995 price from the local Ford Dealership. μ >= 29995
替代假设(H ᴀ ):
Our alternate must be the opposite of the null hypothesis therefore that the average price is not greater then the quoted price from the local Ford dealership. μ < 29995
显著性水平
统计显著性的水平通常表示为所谓的p-值。根据我们稍后选择的统计测试，我们将计算观察样本结果的概率(即p-值)，假设零假设为真。
在进一步讨论之前，让我们先用简单的方法来理解 P 值的含义:
P 值
美国统计协会给出的定义:
在统计假设检验中，p 值或概率值或渐近显著性是给定统计模型的概率，当零假设为真时，统计汇总(如两个比较组之间的样本均值差异)将与实际观察结果相同或比实际观察结果更大。 在学术文献中，p 值被定义为如果零假设为真，数据至少与观察到的一样极端的概率。“至少和观察到的一样极端”——这是什么意思？
为了回答这个问题，让我们用一枚硬币做一个实验。
我的一个朋友给了我一枚硬币，并声称它很特别。我决定测试一下。我的无效假设是:“硬币正常”。然后我扔了 20 次硬币。我有 18 个头和 2 条尾巴。
通常，当扔硬币时，期望是 50%正面和 50%反面，所以我期望 10 正面+ 10 反面的结果是最常见的，其次是 9 正面+ 11 反面和 11 正面+ 9 反面的结果。18 头+ 2 尾的结果到了概率曲线的外围(也就是更极端)。**p 值是观察结果加上所有“更极端”结果的概率(本例中:19 头+ 1 尾和 20 头+ 0 尾)，用阴影“尾面积”**表示。
通常，统计学家会计算一个双尾 p 值(我们稍后会谈到这一点)，所以我们取了图表的两个极端，即 18 个尾+ 2 个头，19 个尾+ 1 个头和 20 个尾+ 0 个头也包括在内。
P 值=数据至少与观察到的数据一样极端的概率= [[ p (18 头 2 尾)+ p (19 头 1 尾)+ p (20 头 0 尾)]+[P(18 尾+ 2 头)+ p (19 尾+ 1 头)+ p (20 尾+ 0 头)]] = 0.0004
** [[]]代表每条尾
获得这样一个结果的机会是如此之小，如果硬币是正常的。因此，我拒绝零假设，并接受硬币是特殊的。
在大多数分析中，0.05 的α被用作显著性的截止值，我们也将在我们的案例中这样做，即如果 p 值小于 0.05，我们将拒绝零假设。
类型 1 &类型 2 错误
理解第一类错误
第一类错误——通常称为“假阳性”——发生在假设检验中，当**零假设为真但被拒绝时。**简单地说，当测试人员验证了统计上的显著差异时，即使没有差异，也会发生这种情况。第一类错误的例子是，当我们推断上面的硬币是特殊的时候，它实际上是正常的，而我们得到 18 个正面和 2 个反面只是运气不好。
类型 1 错误的概率为“α”,与您设置的置信度相关。置信度为 95%的测试意味着有 5%的几率出现 1 类错误。这也是我们面临的风险，5%的机会可能对我们不利。
理解第二类错误
第二类错误被称为“假阴性”，当零假设为假，并且我们随后未能拒绝它时，这些错误就会发生。第二类错误的例子是，如果上面的硬币实际上是特殊的，但我们未能拒绝它是正常的。
单尾或双尾测试
我们常常不得不决定我们的统计检验应该是单尾检验还是双尾检验(也分别称为“定向”和“非定向”检验)。那么，这两者到底有什么区别呢？首先，了解术语“尾巴”在上下文中的含义可能会有所帮助。
尾部是指您正在进行的特定分析的检验统计分布的末端。例如，t 检验使用 t 分布，方差分析使用 f 分布。检验统计量的分布可以有一个或两个尾部，这取决于其形状(见下图)。图中分布的黑色阴影区域是尾部。像 t 和 z 分布这样的对称分布有两个尾部。像 f 分布和卡方分布这样的非对称分布只有一个尾部。这意味着方差分析和卡方检验等分析没有“单尾对双尾”选项，因为它们所基于的分布只有一个尾。
但是如果我们正在进行一个双尾分布的 t 检验呢？现在我们必须决定是单尾还是双尾测试最适合你的研究。
我们使用 0.05 的显著性水平，一个双尾测试分配一半的α用于测试一个方向的统计显著性，另一半α用于测试另一个方向的统计显著性。这意味着. 025 在我们的测试统计分布的每一个尾部。
单尾检验为检测效应提供了更大的能力，因为整个权重只分配给一个方向，正因为如此，每当你对效应的方向有一个假设时，可能会有使用单尾检验的诱惑。但在这样做之前，我们必须考虑错过另一个方向的效果的后果。
当使用双尾检验时，不管假设的关系的方向如何，我们都会检验两个方向上关系的可能性。
例如，我们可能希望使用 t 检验将样本的平均值与给定值 x 进行比较。我们的零假设是平均值等于 x 。**如果平均值明显大于 x 并且平均值明显小于 x ，则双尾测试将同时进行测试。**如果检验统计量位于其概率分布的前 2.5%或后 2.5%，导致 p 值小于 0.05，则认为平均值与 x 显著不同。
但是在我们的例子中，零假设表明我们公司认为 2007 年福特探险者的平均成本大于或等于当地福特经销商的报价 29，995 美元。我们没有兴趣看它是否低于报价，这就是为什么我们将使用单尾检验。
现在我们有了 P 值和单/双尾检验的概念，我们将假设的显著性水平设置为 0.05，作为拒绝真假设的风险，即 1 型错误。我们选择 0.05 作为我们的显著性水平，因为传统上 0.05 用于消费者研究时确定显著性水平。
统计测试
如果你按照下面的图表(查看底部的图表一次)选择测试是很容易的，因为我们不知道整个人口的标准偏差(我们只有性病)。戴夫。对于我们的样本)我们必须使用 t 检验。
当使用一个总体均值的检验统计量时， 在两种情况下，您必须使用 t 分布而不是 Z 分布。 第一种情况是 样本量很小 (低于 30 左右)，第二种情况是当总体标准差、 σ未知 时，你必须使用样本标准差 s 来估计它。在这两种情况下，你都没有更可靠的信息来作为你的结论的基础，所以你必须通过使用 t 分布来为此付出代价
涉及 t 分布的总体均值假设检验称为 t 检验。在这种情况下，检验统计量的公式为:
其中 t ₙ ₋₁是具有 n-1 个自由度的 t 分布的值。
总结我们目前掌握的所有数据:
+===========+====================+
| 29,995.00 | hypothesized value |
+-----------+--------------------+
| 30,474.80 | mean Data |
+-----------+--------------------+
| 1,972.59 | sample std. dev. |
+-----------+--------------------+
| 441.085 | std. error |
+-----------+--------------------+
| 20 | n |
+-----------+--------------------+
| 19 | degrees of freedom |
+-----------+--------------------+
从上式计算 t 值→ 1.087772。p 的对应值为 0.145151。
你可以找到许多在线表格和计算器来计算 p 值和 t 值。
p 值大于α 0.05。在这种情况下，我们无法拒绝零假设。当这种情况发生时，我们说这个结果没有统计学意义。换句话说，我们有理由相信，我们观察到的数据可以单独用偶然来解释。通过解释测试结果，我们证明了无效假设。
测试类型
下面的流程图描述了你可以做哪些类型的测试来确认/拒绝你的零假设。我怎么也找不到这张图表的参考文献，我早就把它保存起来了。请让我知道，如果你能找到参考，我将非常乐意添加它。
根据你的数据类型、抽样方法和假设类型，你必须选择合适的测试。每个测试的过程都与我们在这里展示的相似。如果你想知道上述不同类型测试的细节和用例，有 9 种类型的测试已经在这个链接中描述了，以及它们相应的 R 代码。真的很有帮助。
Flowchart to determine the type of statistical test to be done
参考文献:
- https://statistics . laerd . com/statistical-guides/hypothesis-testing . PHP
- https://dzone.com/articles/what-is-p-value-in-layman-terms
- https://www . students 4 best evidence . net/p-value-in-plain-English-2/
- https://upload . wikimedia . org/Wikipedia/en/0/00/P-value _ graph . png
- https://www . machine learning plus . com/statistics/statistical-significance-tests-r/
- https://www.abtasty.com/blog/type-1-and-type-2-errors/
- https://www . dummies . com/education/math/statistics/how-to-use-the-t-test-to-handle-small-samples-and-unknown-standard-deviations/
- https://www . thoughtco . com/the-difference-between-alpha-and-p-values-3126420
- https://brain mass . com/statistics/hypothesis-testing/hypothesis-testing-real-life-example-117671
- https://www . statistics solutions . com/should-you-use-a-one-tailed-test-or-a-two-tailed-test-for-your-data-analysis/
- https://stats . idre . UCLA . edu/other/mult-pkg/FAQ/general/FAQ-one-tailed-and-two-tailed-tests 之间的区别是什么/
日志——R 中线性和多项式回归实用指南
这是 r 中线性和多项式回归的实用指南。我试图用 King County 数据集涵盖理论和实际实现的基础。
我在数据科学领域的旅程中相对来说还是个新手，这篇文章算是对我所做主题的一个注解。这些可能对其他爱好者有所帮助。
请注意，这些是我的学习笔记，有些部分是从其他来源引用的。我在所有帮助过我的资料的末尾都给出了详细的参考。
回归方程
简单的线性回归精确地估计当 X 改变一定量时 Y 将改变多少。有了相关系数，变量 X 和 Y 可以互换。通过回归，我们尝试使用线性关系(即直线)从 X 预测 Y 变量:
Y = b0 + b1 X
符号 b0 称为截距(或常数)，b1 表示 x 的斜率，二者都以系数形式出现在 R 输出中。Y 变量被称为响应变量或因变量，因为它依赖于 X。X 变量被称为预测变量或自变量。机器学习社区倾向于使用其他术语，称 Y 为目标，X 为特征向量。
拟合值和残差
回归分析中的重要概念是拟合值和残差。一般来说，数据不会正好落在一条线上，因此回归方程应该包括一个明确的误差项:
拟合值也称为预测值，通常用(Y-hat)表示。这些由以下公式给出:
其中 b0 & b1 表示系数是估计的还是已知的。
多元线性回归
当有多个预测值时，该等式可以简单地扩展以适应它们:
我们现在有了一个线性模型，而不是直线，即每个系数与其变量(特征)之间的关系是线性的。我们将有一个平面，它将通过“中间”的采样点，如下图所示。二维是一条线，三维是一个平面，N 是一个超平面。
现在我们对什么是线性回归有了一个概念，我们将通过不同类型回归模型的细节来预测非常受欢迎的 King County 住房数据集的价格。您可以从 这里 下载数据集。目标是我们使用不同的模型来定义一个平面，该平面以低错误率最适合我们的数据。随着新数据的出现，该平面/线将用于预测结果。
郡王房产数据
在开始分析之前，请确保您已经安装并加载了以下 R 包:
library(doParallel)
library(‘randomForest’)
library(MASS)
library(readr)
library(ggplot2)
library(corrplot)
library(mlbench)
library(Amelia)
library(plotly)
library(reshape2)
library(caret)
library(caTools)
library(dplyr)
library(moments)
探索性数据分析
使用回归的一个例子是估计房屋的价值。郡估税员必须估计房子的价值，以达到估税的目的。房地产消费者和专业人士会咨询 Zillow 等热门网站，以确定一个合理的价格。以下是来自华盛顿州金县(西雅图)的几行住房数据，来自住房数据框:
house <- read.csv(“C:/personal/kc_house_data.csv”)
head(house)
数据集本身由 21，613 个具有 19 个特征的示例组成[不考虑 id 和日期列]。其特点是:
现在，如果你注意到这些数据，你会发现很少有值是被分类的，比如海滨、风景、条件、等级等等。现在还有一个对邮政编码变量进行分组的空间，因为邮政编码本身并不定义任何东西，但是房子的价格取决于位置。如果你知道一些邮政编码属于城镇的豪华区，那么你可以把它们分成邮政编码组 1，下一个是邮政编码组 2，依此类推。我们不会在这里这样做，因为我们不知道哪个邮政编码区属于城镇的哪个部分。但是我们有 lat 和 long 值，所以我们将使用它们。但是如果你想了解更多关于将数据转换成分类值的信息，你可以参考这个链接。通常使用电平命令来完成。
缺失值
您可以使用以下命令检查丢失的数据:
missmap(house,col=c(‘yellow’,’black’),y.at=1,y.labels=’’,legend=TRUE)
根据该图，没有缺失值。
现在，我们将删除 id 列和日期列(所有信息仅在一年内收集)。
house <- house[,-c(1,2)]
相关预测因子
在多元回归中，预测变量通常是相互关联的。相关性是一个问题，因为独立变量应该是独立的，如果两个变量相关，那么这意味着你基本上多次使用几乎相同的特性。这将给出不正确的结果。一旦你找到了相关变量，你就可以保留一个，去掉另一个。
如果删除高度相关的属性，许多方法会执行得更好。您可以使用以下命令检查相关变量:
corrplot(cor(house),type=”upper”,method=”color”,addCoef.col = “black”,tl.col = “black”, number.cex = 0.6,mar=c(0,0,0,0))
现在我们知道了哪些变量对是相关的，我们可以选择保留哪些变量对，删除哪些变量对:
[sqft_living, grade] : 0.76 , will keep sqft_living
[sqft_living, sqft_above] : 0.88 , will keep sqft_living
[sqft_living, sqft_living15] : 0.76 , will keep sqft_living
[sqft_lot, sqft_lot15] : 0.72 , will keep sqft_lot
因此，我们删除了 grade、sqft_above、sqft_living15 和 sqft_lot15，因为它们所携带的信息已经由其他变量提供了。
house <- house[,-c(10,11,18,19)]
维度
如果您想减少数据集的维度，您可以应用脱字符包中的函数 nearZeroVar 。它诊断相对于样本数具有一个或很少唯一值的预测值，并且最常见值的频率与第二常见值的频率之比很大。截至目前，我们不会这样做，因为我们还有 14 个功能。
虽然我会从数据集中删除邮政编码，因为位置数据包含在纬度和经度中。
house <- house[,-c(13)]
偏斜度
接下来，我们将尝试分析数字变量的偏斜度。进行偏斜度检查是为了使数据正常，这是线性回归的假设之一。有些人建议这里偏斜度的可接受值范围在(-2，2)之间。因此，我们检测哪些变量不在此范围内，并将使用对数函数对它们进行转换[我们不会转换分类数据]。
apply(house, 2, skewness, na.rm =TRUE)
将所需数据转换为日志:
house$price <- log(house$price)
house$sqft_lot <- log(house$sqft_lot)
可以在下面两个直方图中检查转换前后偏斜度的差异:
数据建模
现在我们已经完成了 EDA，我们将从从其他变量预测销售价格的目标开始。首先，我们将采用所有变量/特征来制作我们的模型；参数 na.action=na.omit 导致模型删除缺少值的记录。
使用 70/30 分割将加载的数据集分成训练集和验证集。现在，我们有了一个用于训练模型的训练数据集和一个稍后用于测量模型性能的验证集。
set.seed(1000)
i <- sample(2, nrow(house), replace=TRUE, prob=c(0.7, 0.3))
trainData <- house[i==1,]
testData <- house[i==2,]
model_lm <- lm(price ~ bedrooms + bathrooms + sqft_living + sqft_lot + floors +waterfront+ view + sqft_basement + condition + yr_built + lat + long, data = trainData)
summary(model_lm)
正如我们可以看到运行我们的线性模型的总结，决定系数(或 R 平方)是好的。从我们创建的线性回归模型中可以预测结果中 70.8%的方差。
接下来，我们将分析我们模型的 RMSE(或均方根误差)。
prediction_lm <- predict(model_lm, testData, type=”response”)
model_lm_output <- cbind(testData, prediction_lm)
# RMSE value
RMSE_lm <- sqrt(mean((model_lm_output$price — model_lm_output$prediction_lm)^2,na.rm=TRUE))
print(RMSE_lm)
正如我们在上面看到的，我们的 RMSE 是 0.28。它衡量模型预测的价格和实际价格之间的差异。数值越低越好。我们的接近 0，所以这是一个很好的指标。
现在我们有了一个基本模型，我们需要了解模型的性能:
从数据科学的角度来看，最重要的性能指标是均方根误差，即 RMSE。 RMSE 是预测值的均方误差的平方根:
您可以通过以下方式计算上述模型中的 RMSE
RMSE_house_lm <- sqrt(mean(residuals(house_lm)^2))
这衡量了模型的整体准确性，并且是将其与其他模型(包括使用机器学习技术拟合的模型)进行比较的基础。
与 RMSE 相似的是剩余标准误差，或称 RSE 。在这种情况下，我们有 p 个预测值，RSE 由下式给出:
MSE_house_lm <- (mean(residuals(house_lm)^2))
唯一的区别是分母是自由度，而不是记录数。在实践中，对于线性回归，RMSE 和 RSE 之间的差异非常小，特别是对于大数据应用。
在我们的例子中，观察到的价格值平均偏离预测价格约 202700 个单位。这对应的误差率 202700/mean(房价)= 37.5% ，偏高。
残差平方和 (RSS)是残差平方和:
rss <- sum(residuals(fit)^2)
残差标准差 (RSE)是(RSS /自由度)的平方根:
rse <- sqrt( sum(residuals(fit)^2) / fit$df.residual )
您将在软件输出中看到的另一个有用的度量是决定系数，也称为 R 平方统计量。R-squared 的范围从 0 到 1，用于测量模型中数据的变化比例。它主要用于回归的解释性应用，在这种情况下，您需要评估模型与数据的拟合程度。公式是:
分母与 y 的方差成比例，R 的输出也报告调整后的 R 平方，它针对自由度进行调整；这在多元回归中很少是显著不同的。
在我们的例子中，R2 为 0.69，这是可以的。
除了估计的系数，R 还报告了系数的标准误差(SE)和 t 统计量:
t 统计量(及其镜像 p 值)测量系数的“统计显著性”程度，即超出预测值和目标变量随机排列可能产生的范围。t 统计值越高(p 值越低)，预测值越显著。
警告
除了 t 统计量之外，R 和其他软件包通常会报告 p 值(R 输出中的 Pr(>|t|)和 F 统计量。数据科学家通常不会过多参与这些统计数据的解释，也不会过多参与统计显著性的问题。数据科学家主要关注 t 统计，将其作为是否在模型中包含预测因子的有用指南。较高的 t 统计值(p 值接近 0)表示模型中应保留某个预测值，而非常低的 t 统计值则表示可能会删除某个预测值。
我们可以从线性模型的图示中得到一些启示:
# Graphic representation
par(mfrow=c(2,2))
plot(model_lm)
现在我们如何阅读这些图表？
- 残差与拟合值:由于残差随机分布在零附近，因此不应有强模式和异常值。正如我们在图中看到的，我们所有的残差都分布在(-1，1)内，接近于 0。
- 正态 Q-Q :残差应该正态分布在对角线周围，就像我们的例子中发生的那样。
- 缩放位置:与第一个图类似，但是 Y 轴上的残差被缩放。不应该有任何强烈的模式，因为它发生在这种情况下，虽然我们有一个平滑的正对角线
- 残差 vs 杠杆:检测异常值很有用。在我们的例子中，似乎我们的图不包含任何超出不连续红线(库克距离)的异常值。看看 Cook 的距离值:任何大于 1 的值都表示可能影响模型的情况。在统计学中， Cook’s distance 或 Cook’s D 是执行最小二乘回归分析时对数据点影响的常用估计。在实际的普通最小二乘分析中，库克距离可以以几种方式使用:指示特别值得检查有效性的有影响的数据点；或者指示设计空间中能够获得更多数据点的区域。
一些教科书告诉你，库克距离大于 1 的点被认为是有影响的。其他文本给你的阈值是 4/ N 或 4/(N—k—1)，其中 N 是观察值的个数， k 是解释变量的个数。约翰·福克斯在他关于回归诊断的小册子中，在给出数字阈值时相当谨慎。他建议使用图形，并更仔细地检查“D 值远大于其余值”的点。根据 Fox 的说法，阈值应该只是用来增强图形显示。在我们的数据中，我们可以看到至少有两点值得注意。
cooksd <- cooks.distance(model_lm)
# Plot the Cook’s Distance
sample_size <- nrow(trainData)
plot(cooksd, pch=”*”, cex=2, main=”Influential Obs by Cooks distance”) # plot cook’s distance
text(x=1:length(cooksd)+1, y=cooksd, labels=ifelse(cooksd>0.1, names(cooksd),””), col=”red”) # add labels
使用上面的命令和图形，我们可以获得异常数据的行号，如果您手动检查这些数据，您会看到这些数据中存在差异(一个人有 33 间卧室)，让我们看看如果我们删除这些会发生什么:
#try to remove the top x outliers to have a look
top_x_outlier <- 2
influential <- as.numeric(names(sort(cooksd, decreasing = TRUE)[1:top_x_outlier]))
house2 <- house[-c(15871,12778), ]
再次运行该模型，我们得到了下面的 RMSE，它比之前的结果略有改进:
按重要性排列特征
很难猜测我们的模型中包括哪些特征，因为我们不知道选择哪些特征，我们在初始模型中采用了所有特征。通过建立模型，可以从数据中估计特征的重要性。像决策树这样的方法有一个内置的机制来报告可变的重要性。对于其他算法，可以使用对每个属性进行的 ROC 曲线分析来估计重要性。然后，varImp 用于估计变量重要性，并打印和绘制出来。
# prepare training scheme
control <- trainControl(method=”repeatedcv”, number=10, repeats=3)
# train the model
model <- train(price~., data=trainData, method=”lm”, preProcess=”scale”, trControl=control)
# estimate variable importance
importance <- varImp(model, scale=FALSE)
# summarize importance
print(importance)
# plot importance
plot(importance)
运行上面的命令后，我们可以看到模型认为重要的特性。
功能选择
自动特征选择方法可用于利用数据集的不同子集构建许多模型，并识别构建精确模型所需的和不需要的那些属性。
由 caret R 包提供的一种流行的自动特征选择方法被称为递归特征消除或 RFE。
以下示例提供了数据集上 RFE 方法的示例。在每次迭代中使用随机森林算法来评估模型。该算法被配置成探索所有可能的属性子集。本例中选择了所有 13 个属性，尽管在显示不同属性子集大小的图中，我们可以看到只有 6 个属性给出了几乎相当的结果。
现在，RFE 算法的问题是运行需要大量的时间(我说的是 10-12 个小时)，为了减少时间，我们可以利用集群(减少到近 2 个小时，在此期间尝试检查您的任务管理器)。
要使用多个工作线程调整预测模型，需要使用一个单独的函数来“注册”并行处理技术，并指定要使用的工作线程数量。例如，要在同一台机器上使用带有五个内核的 doParallel 包，需要加载并注册该包(检查一次任务管理器，您会注意到不同之处):
library(doParallel)
cl <- makePSOCKcluster(5)
registerDoParallel(cl)
control <- rfeControl(functions=rfFuncs, method=”cv”, number=10)
# run the RFE algorithm, once the cluster cl is built then the algorithm is by default run in parallel.
results <- rfe(house[,2:19], house[,1], sizes=c(1:19), rfeControl=control)
# summarize the results
print(results)
# list the chosen features
predictors(results)
# plot the results
plot(results, type=c(“g”, “o”))
# once the job is done, stop the cluster
stopCluster(cl)
根据上述模型输出，以下 6 个变量产生最低的 RMSE:
“lat” “long” “sqft_living” “sqft_lot” “yr_built” “view”
如果您注意到在 RFE 和 varIMP 方法中选择的特征有一点不同，不同的方法将选择不同的特征子集。就像没有最好的模型一样，可能没有“最好的”特性集。
交叉验证
现在，如果你已经注意到我们正在使用 method=‘CV’ ，那么接下来我们将看看这个 CV 代表什么。经典的统计回归度量(R2、F 统计和 p 值)都是“样本”度量，它们应用于用于拟合模型的相同数据。直观地说，您可以看到，将一些原始数据搁置起来，而不是用它来拟合模型，然后将模型应用于搁置(维持)数据，以查看它的效果，这是很有意义的。通常，您会使用大部分数据来拟合模型，并使用一小部分数据来测试模型。这种“样本外”验证的想法并不新鲜，但直到更大的数据集变得更加普遍，它才真正站稳脚跟；对于小型数据集，分析师通常希望使用所有数据并拟合最佳模型。但是，使用维持样本会使您受到一些不确定性的影响，这些不确定性仅仅是由小型维持样本的可变性引起的。如果您选择不同的维持样本，评估会有多大的不同？
交叉验证将维持样本的概念扩展到多个连续的维持样本。基本 k 倍交叉验证的算法如下:
1\. Set aside 1/k of the data as a holdout sample.
2\. Train the model on the remaining data.
3\. Apply (score) the model to the 1/k holdout, and record needed model assessment metrics.
4\. Restore the first 1/k of the data, and set aside the next 1/k (excluding any records that got picked the first time).
5\. Repeat steps 2 and 3.
6\. Repeat until each record has been used in the holdout portion.
7\. Average or otherwise combine the model assessment metrics.
将数据分为定型样本和维持样本也称为折叠。
逐步回归
在我们的问题中，许多变量可以在回归中用作预测因子。有一种称为逐步回归的方法，它可以从模型中依次添加或删除特征，并检查模型的性能。
请注意，train()函数[caret package]提供了一个简单的工作流来使用 leaps 和 MASS 包执行逐步选择。它有一个名为 method 的选项，可以采用以下值:
- “leapBackward”，用向后选择拟合线性回归
- “前向”，用前向选择拟合线性回归
- “leapSeq”，用逐步选择拟合线性回归。
您还需要指定调整参数 nvmax，它对应于要合并到模型中的最大预测值数量。
例如，您可以在 1 到 5 之间改变 nvmax。在这种情况下，该函数从搜索不同大小的不同最佳模型开始，直到最佳 5 变量模型。也就是说，它搜索最佳单变量模型，最佳双变量模型，…，最佳五变量模型。
我们将使用 10 重交叉验证来估计 5 个模型的平均预测误差(RMSE)。RMSE 统计指标用于比较 5 个模型，并自动选择最佳模型，其中最佳定义为最小化 RMSE 的模型。
逐步回归对于包含多个预测变量的高维数据非常有用。其他替代方法是惩罚回归(岭和套索回归)和基于主成分的回归方法(PCR 和 PLS)。
## THIS IS FOR STEPWISE REGRESSION
set.seed(123)
# Set up repeated k-fold cross-validation
train.control <- trainControl(method = “cv”, number = 10)
# Train the model
step.model <- train(price ~., data = trainData[,1:14],method = “leapSeq”,tuneGrid = data.frame(nvmax = 1:13),trControl = train.control)step.model$results
step.model$bestTune
summary(step.model$finalModel)
有 7 个变量的 RMSE 几乎和我们一开始发现的有更多变量的那个相似。
增加更多的变量，并不一定意味着我们有一个更好的模型。统计学家使用奥卡姆剃刀原理来指导模型的选择:在所有条件相同的情况下，应优先使用更简单的模型，而不是更复杂的模型。
包括额外的变量通常会减少 RMSE，增加 R2。因此，这些不适合帮助指导模型选择。20 世纪 70 年代，日本杰出的统计学家 Hirotugu Akaike 开发了一种叫做 AIC (Akaike 的信息标准)的指标，对在模型中添加术语进行惩罚。在回归的情况下，AIC 具有以下形式:
AIC = 2P + n log(RSS/n)
其中 p 是变量的数量，n 是记录的数量。目标是找到最小化 AIC 的模型；多 k 个额外变量的模型被罚 2k。
AIC、BIC 和锦葵 CP
AIC 的公式似乎有点神秘，但事实上它是基于信息论的渐近结果。AIC 有几种变体:
AICc:针对小样本量进行修正的 AIC 版本。
BIC 或贝叶斯信息标准:类似于 AIC，但在模型中加入了额外的变量，惩罚力度更大。
锦葵 Cp:由科林·锦葵开发的 AIC 的变体。
数据科学家通常不需要担心这些样本内指标之间的差异或它们背后的基本理论。
我们如何找到最小化 AIC 的模型？一种方法是搜索所有可能的模型，称为所有子集回归。这在计算上是昂贵的，并且对于具有大量数据和许多变量的问题是不可行的。一个有吸引力的替代方法是使用逐步回归，我们在上面学过，这种方法连续增加和减少预测因子来找到降低 AIC 的模型。更简单的是向前选择和向后选择。在正向选择中，您从没有预测值开始，然后逐个添加它们，在每一步添加对贡献最大的预测值，当贡献不再具有统计显著性时停止。在向后选择或向后排除中，从完整的模型开始，去掉不具有统计显著性的预测因子，直到剩下一个所有预测因子都具有统计显著性的模型。
逐步方法的问题在于，它们基于模型中的其他变量来评估变量的拟合度。有人用穿衣服的比喻来形容这个问题。如果一个逐步回归方法正在选择你的衣服，它会根据已经选择的衣服决定你应该穿什么衣服。
惩罚退化在精神上类似于 AIC。模型拟合方程结合了一个约束，该约束由于太多的变量(参数)而对模型不利，而不是显式地搜索一组离散的模型。惩罚回归不是像逐步选择、向前选择和向后选择那样完全消除预测变量，而是通过降低系数(在某些情况下接近于零)来施加惩罚。常见的惩罚回归方法有岭回归和拉索回归。
逐步回归和所有子集回归是评估和调整模型的样本内方法。这意味着模型选择可能会过度拟合，并且在应用于新数据时可能表现不佳。避免这种情况的一种常见方法是使用交叉验证来验证模型。在线性回归中，过度拟合通常不是主要问题，因为对数据施加了简单的(线性)全局结构。对于更复杂类型的模型，特别是响应局部数据结构的迭代过程，交叉验证是一个非常重要的工具。
山脊和套索回归
本节摘自本文 优析 vidhya 文章， 要了解更多关于脊和套索回归背后的数学请做浏览链接。
正则化有助于解决机器学习中的过拟合问题。简单的模型将是一个非常糟糕的数据概括。同时，由于过度拟合，复杂模型在测试数据中可能表现不佳。我们需要在简单模型和复杂模型之间选择合适的模型。正则化有助于选择首选的模型复杂性，以便模型更好地预测。正则化无非是在目标函数中增加一个惩罚项，并使用该惩罚项来控制模型的复杂性。它可以用于许多机器学习算法。岭和套索回归都使用 L2 和 L1 正则化。
岭和套索回归是强大的技术，通常用于在存在“大量”要素的情况下创建简约的模型。这里的“大”通常意味着两种情况之一:
- 大到足以增强模型过度拟合的趋势(低至 10 个变量就可能导致过度拟合)
- 大到足以导致计算挑战。对于现代系统，这种情况可能会出现在数百万或数十亿个特征的情况下
尽管山脊线和套索看起来似乎是为了一个共同的目标而工作，但它们的内在属性和实际使用情况却大相径庭。如果您以前听说过它们，您一定知道它们通过惩罚特征系数的大小以及最小化预测和实际观测之间的误差来工作。这些被称为“正则化”技术。关键区别在于他们如何给系数分配惩罚:
- 岭回归:
- 执行 L2 正则化，即添加相当于系数大小的平方的惩罚
- 最小化目标= LS Obj + α *(系数的平方和)
2。拉索回归:
- 执行 L1 正则化，即添加相当于系数大小的绝对值的惩罚
- 最小化目标= LS Obj + α *(系数绝对值之和)
注意，这里的“LS Obj”指的是“最小二乘目标”，即没有正则化的线性回归目标。
现在我们对岭和套索回归的工作原理有了一点了解，让我们试着通过比较它们来巩固我们的理解，并试着欣赏它们的具体用例。我还会将它们与一些替代方法进行比较。让我们从三个方面来分析这些问题:
1。关键区别
- **脊:**它包括模型中的所有(或没有)特征。因此，岭回归的主要优点是系数收缩和降低模型复杂性。
- **套索:**除了收缩系数，套索还执行特征选择。一些系数变得恰好为零，这相当于从模型中排除了特定的特征。
传统上，像逐步回归这样的技术被用来执行特征选择和制作简约的模型。但是随着机器学习的进步，脊和套索回归提供了非常好的选择，因为它们提供了更好的输出，需要更少的调整参数，并且可以在很大程度上自动化。
2。典型使用案例
- **脊:**主要用于防止过配合。因为它包括所有的特征，所以它在过高的特征的情况下不是很有用，比如以百万计，因为它会带来计算上的挑战。
- Lasso: 因为它提供了稀疏解，所以它通常是#特征数以百万计或更多的建模案例的首选模型(或该概念的某种变体)。在这种情况下，获得稀疏解具有很大的计算优势，因为可以简单地忽略具有零系数的特征。
不难看出为什么逐步选择技术在高维情况下变得非常难以实施。因此，套索提供了一个显著的优势。
3。存在高度相关的特征
- **岭:**即使存在高度相关的特征，它通常也能很好地工作，因为它会将所有特征包括在模型中，但系数将根据相关性在它们之间分配。
- Lasso: 从高度相关的特征中任意选择一个特征，将其余特征的系数降为零。此外，所选择的变量随着模型参数的变化而随机变化。与岭回归相比，这种方法通常不太有效。
除了脊和套索，弹性网格是另一种结合了 L1 和 L2 正则化的有用技术。它可以用来平衡脊和套索回归的利弊。我鼓励你进一步探索它。
现在让我们在数据集上尝试岭套索回归:
#Ridge Regression Modeltr.control <- trainControl(method=”repeatedcv”, number = 10,repeats = 10)
lambdas <- seq(1,0,-.001) #Tuning Parameter
set.seed(123)
ridge_model <- train(price~., data=trainData,method=”glmnet”,metric=”RMSE”,maximize=FALSE,trControl=tr.control,tuneGrid=expand.grid(alpha=0,lambda=lambdas))
ridge_model$results
varImp(ridge_model)
ridge_preds <- predict(ridge_model,newdata = testData)
model_rrm_output <- cbind(testData, ridge_preds)
RMSE_rrm <- sqrt(mean((model_rrm_output$price - model_rrm_output$ridge_preds)^2,na.rm=TRUE))print(RMSE_rrm)
接下来，我们将尝试套索回归，看看我们是否能从中获得更好的输出。
#Lasso Regression Modelset.seed(123)
lasso_model <- train(price~., data=trainData,method=”glmnet”,metric=”RMSE”,maximize=FALSE,trControl=tr.control,tuneGrid=expand.grid(alpha=1,lambda=c(1,0.1,0.05,0.01,seq(0.009,0.001,-0.001), 0.00075,0.0005,0.0001)))varImp(lasso_model)
lasso_preds <- predict(lasso_model,newdata = testData)
model_lrm_output <- cbind(testData, lasso_preds)RMSE_lrm <- sqrt(mean((model_lrm_output$price - model_lrm_output$ridge_preds)²,na.rm=TRUE))print(RMSE_lrm)
请记住 RMSE 分数，我们将使用它来进一步比较模型性能。
加权回归
统计学家出于各种目的使用加权回归；特别是，它对于复杂调查的分析非常重要。数据科学家可能会发现加权回归在两种情况下很有用:
- 当不同的观测值以不同的精度测量时，采用逆方差加权。
- 以聚集形式对数据进行分析，使得权重变量对聚集数据中每行代表的原始观察值进行编码。
例如，对于住房数据，旧的销售数据不如最近的销售数据可靠。
既然我们已经尝试了不同的线性回归模式，我们可以转向其他方法，看看我们是否可以进一步降低 RMSE。
非线性回归
在某些情况下，结果和预测变量之间的真实关系可能不是线性的。
有不同的解决方案来扩展线性回归模型以捕捉这些非线性效应，包括:
- 多项式回归。这是建模非线性关系的简单方法。它将多项式项或二次项(平方、立方等)添加到回归中。
- 样条回归。用一系列多项式段拟合平滑曲线。界定样条线段的值称为节点。
- 广义可加模型 (GAM)。通过自动选择结来拟合样条模型。
假设您怀疑响应和预测变量之间的非线性关系，无论是通过先验知识还是通过检查回归诊断。多项式项可能不够灵活，无法捕捉这种关系，而样条项需要指定节点。广义加性模型(GAM)是一种自动拟合样条回归的技术。因此，我们接下来将在我们的数据上尝试 GAM。
广义可加模型
如果你绘制价格与纬度的关系图，你会发现多项式曲线比直线拟合得更好，你可以对其他变量进行同样的分析，并检查它们的行为。
plot(house$lat,house$price)
下面的回归曲线拟合得更好:
thegam<-gam(price ~ s(lat),data=trainData)
plot(thegam, residuals=TRUE,shade.col=”gray80")
对我们的数据运行 GAM，我们有:
library(mgcv)lm_gam <- gam(price ~ bedrooms + bathrooms + s(sqft_living) + sqft_lot + floors +waterfront+ view + s(sqft_basement) + condition + yr_built + s(lat) + s(long),data=trainData)prediction_lm <- predict(lm_gam, testData, type=”response”)model_lm_output <- cbind(testData, prediction_lm)# RMSE value
RMSE_lm <- sqrt(mean((model_lm_output$price-model_lm_output$prediction_lm)²,na.rm=TRUE))print(RMSE_lm)
术语 s(sqfttoliving)告诉 gam 函数找到样条术语的“最佳”节点。
你可以看到 RMSE 分数比正常回归要好。
随机森林
随机森林是一种能够执行回归和分类任务的算法。在回归的情况下，它通过在训练时构建大量决策树并输出作为各个树的平均预测的类来操作。我们将在后面详细研究随机森林算法。
正如我们之前在线性模型中所做的那样，我们将数据集分为训练集和验证集。之后，我们定义模型中包含的变量并运行它。请注意，我们将采用该模型中的所有变量:
library(randomForest)set.seed(2000)tree_rf <- price ~ bedrooms + bathrooms + sqft_living + sqft_lot + floors + condition + grade + sqft_above + yr_built + zipcode + lat + long + sqft_living15 + sqft_lot15model_rf <- randomForest(tree_rf, data=trainData, ntree=100,proximity=T, do.trace=10)varImpPlot(model_rf)
plot(model_rf, main = “”)# Making predictions
prediction_rf <- predict(model_rf, testData)
model_rf_output <- cbind(testData, prediction_rf)# RMSE
RMSE_rf <- sqrt(mean((model_rf_output$price - model_rf_output$prediction_rf)²,na.rm=TRUE))
正如我们在上面看到的，随机森林模型的 RMSE 是 0.18。RMSE 值是迄今为止我们看到的最低值。
渐变助推
梯度推进是解决机器学习问题的最有力的技术之一。这是一种集成学习算法，它结合了几个基本估计量的预测，以提高单个估计量的鲁棒性。作为随机森林，它也是一个树分类器模型。稍后我们将了解更多关于 GBM 的内容。
以下代码对我们的数据运行 GBM:
fitControl <- trainControl(method = “cv”, number = 50)tune_Grid <- expand.grid(interaction.depth = 2, n.trees = 500, shrinkage = 0.1, n.minobsinnode = 10)set.seed(825)tree_gbm <- price ~ bedrooms + bathrooms + sqft_living + sqft_lot + floors + condition + grade + sqft_above + yr_built + zipcode + lat + long + sqft_living15 + sqft_lot15model_gbm <- train(tree_gbm, data = trainData, method = “gbm”, trControl = fitControl, verbose = FALSE)print(model_gbm)
prediction_lm <- predict(model_gbm, testData, type=”response”)
model_lm_output <- cbind(testData, prediction_lm)# RMSE value
RMSE_lm <- sqrt(mean((model_lm_output$price-model_lm_output$prediction_lm)²,na.rm=TRUE))print(RMSE_lm)
这也不错，但不如随机森林结果好。总结到目前为止我们在本文中看到的所有模型的结果:
结果很明显随机森林赢了。但是请注意，RF 并不一定每次都表现得更好，对于不同的数据集，不同的模型会优于 RF，我们的工作是找到、调整和使用提供最佳结果的模型。
这篇文章很长，但希望你能从中有所收获。请仔细阅读下面的参考资料，它们包含了我在这里无法捕捉到的主题的更多细节。下次见。
引用:
- https://ccapella . github . io/post/predicting-king-county-USA-house-prices/
- https://www . ka ggle . com/anwarmohammad/housing-prices-prediction-using-r-an war/code
- https://stats . stack exchange . com/questions/110999/r-confused-on-residual-terminals
- https://stats . stack exchange . com/questions/22161/how-to-read-cooks-distance-plots
- https://topepo.github.io/caret/parallel-processing.html
- http://www . sth da . com/English/articles/37-model-selection-essentials-in-r/154-stepwise-regression-essentials-in-r/
- https://courses . analyticsvidhya . com/courses/introduction-to-data-science-2？UTM _ source = blog&UTM _ medium = rideandlassegregationarticle
- http://www . sth da . com/English/articles/40-regression-analysis/162-nonlinear-regression-essentials-in-r-polynomial-and-spline-regression-models/
- https://www . analyticsvidhya . com/blog/2016/01/complete-tutorial-ridge-lasso-regression-python/
书籍:
- 统计学习及其应用导论
- 100 页的机器学习书籍——安德烈·布尔科夫
- 使用 r-Andy Field | Jeremy Miles | zo Field 发现统计数据
- 数据科学家实用统计学— Peter Bruce | Andrew Bruce
- 实用数据科学与 R-NINA ZUMEL | JOHN MOUNT
Apache Kafka 中的日志压缩主题
当我开始阅读 Kafka 文档时，尽管日志压缩主题似乎是一个简单的概念，但我并不清楚 Kafka 如何在文件系统中内部保存它们的状态。这个月我有时间阅读了更多关于这个特性的内容，我想和你分享我的理解。
TL；速度三角形定位法(dead reckoning)
在这篇文章中，我将描述卡夫卡中的日志压缩主题。然后，我将向您展示 Kafka 如何在文件系统中内部保存这些主题的状态。
先决条件
我假设您已经熟悉 Apache Kafka 的基本概念，如代理、主题、分区、消费者和生产者。另外，如果您想运行示例命令，您必须运行一个 Kafka broker 和一个 Zookeeper 服务器。
什么是日志压缩主题
卡夫卡的文献记载说:
日志压缩是一种机制，可以提供更细粒度的按记录保留，而不是更粗粒度的基于时间的保留。这个想法是有选择地删除具有相同主键的更新记录。这样可以保证日志中至少有每个键的最新状态。
为了简化这种描述，当分区日志中存在具有相同关键字的新版本时，Kafka 会删除任何旧记录。作为示例，考虑名为“最新产品价格”的日志压缩主题的以下部分:
正如您首先看到的，有两条带有键 p3 的记录。但是因为这是一个日志压缩的主题，Kafka 在一个后台线程中删除了旧的记录(在接下来的小节中会有更多的介绍)。现在假设我们有一个生产者向这个分区发送新记录。生产者产生 3 个记录，分别具有键 p6、p5、p5:
Kafka broker 中的一个后台线程再次删除了带有键 p5 和 p6 的旧记录。请注意，压缩日志由两部分组成:尾部和头部。Kafka 确保尾部分中的所有记录都有唯一的键，因为尾部分是在清理过程的前一个循环中扫描的。但是 head 部分可以有重复的值。
既然我们已经了解了什么是日志压缩主题，那么是时候使用kafka-topics工具来创建它们了。
创建日志压缩主题
创建一个压缩主题(我将详细描述所有配置):
kafka-topics --create --zookeeper zookeeper:2181 --topic latest-product-price --replication-factor 1 --partitions 1 --config "cleanup.policy=compact" --config "delete.retention.ms=100" --config "segment.ms=100" --config "min.cleanable.dirty.ratio=0.01"
制作一些记录:
kafka-console-producer --broker-list localhost:9092 --topic latest-product-price --property parse.key=true --property key.separator=:
>p3:10$
>p5:7$
>p3:11$
>p6:25$
>p6:12$
>p5:14$
>p5:17$
注意，在上面的命令中，我用:分隔了键和值。现在消费题目:
kafka-console-consumer --bootstrap-server localhost:9092 --topic latest-product-price --property print.key=true --property key.separator=: --from-beginning
p3:11$
p6:12$
p5:14$
p5:17$
如您所见，带有重复键的记录被删除。p5:14$记录没有被删除，我们将在描述清理过程时看到原因。但我们首先要看卡夫卡内部是如何储存信息的。
片段
分区日志是一种抽象，它允许我们轻松地消费分区内部的有序消息，而不必担心 Kafka 的内部存储。然而实际上，分区日志被 Kafka broker 分成了段。段是存储在文件系统中的文件(在数据目录和分区目录中)，它们的名称以.log结尾。在下图中，分区日志分为 3 个部分:
如您所见，我们有一个分区日志，它包含 7 条记录，分别位于 3 个独立的段文件中。一个段的第一个偏置称为该段的基准偏置。段文件名始终等于其基本偏移值。
分区中的最后一个段称为活动段。只有日志的活动段可以接收新生成的消息。我们将会看到卡夫卡在清理压缩原木的过程中是如何处理活动段的。
回到我们的例子，我们可以通过下面的命令查看我们的主题分区的段文件(假设您的 Kafka 数据目录是/var/lib/kafka/data):
ls /var/lib/kafka/data/latest-product-price-0/
00000000000000000000.index 00000000000000000006.log
00000000000000000000.log 00000000000000000006.snapshot
00000000000000000000.timeindex 00000000000000000006.timeindex
00000000000000000005.snapshot leader-epoch-checkpoint
00000000000000000006.index
00000000000000000000.log和00000000000000000006.log是该分区的段，00000000000000000006.log是活动段。
卡夫卡什么时候创造了一个新的片段？一种方法是在主题创建期间设置segment.bytes(默认为 1GB)配置。当您的分段大小变得大于该值时，Kafka 将创建一个新的分段。另一个选择是通过设置segment.ms，就像你之前看到的那样。使用该选项，当 Kafka 收到一个生产请求时，它将检查活动段是否比segment.ms值旧。如果它是旧的，那么它将创建一个新的段。在我们的命令中，我们设置segment.ms=100 来确保每 100 毫秒创建一个新的段。
有趣的是，当你设置segment.ms=100时，你可能会有更小的段。在清理过程(见下一节)之后，Kafka broker 将合并非活动段并从中创建一个大段。
有关 Kafka 的分段和内部存储的更多信息，您可以阅读Kafka 的存储内部结构如何工作和Kafka 存储内部结构实用介绍文章。
洁化过程
在启动过程中，Kafka broker 会创建一些清理线程，负责清理压缩的日志(这些线程的数量可以通过log.cleaner.threads config 进行配置)。清洁线程会不断尝试在代理中找到最脏的日志，然后尝试清理它。对于每个日志，它按如下方式计算脏比:
dirty ratio = the number of bytes in the head / total number of bytes in the log(tail + head)
cleaner 线程然后选择具有最高脏比的日志。这个日志被称为最脏的日志，如果它的值大于min.cleanable.dirty.ratio config，它将被清除。否则，清洁线程将被阻塞数毫秒(可通过log.cleaner.backoff.ms配置)。
在找到最脏的日志后，我们希望找到可清理的日志部分。请注意，日志的某些部分是不可清除的，因此不会被扫描:
- 活动段内的所有记录。这就是为什么我们仍然在消费者中看到重复的 p5:14$记录。
- 如果您将
min.compaction.lag.ms配置设置为大于 0，那么任何包含时间戳早于此配置的记录的段都不会被清除。不会对这些段进行压缩扫描。
现在我们知道了要压缩哪些记录。从日志中的第一条记录到第一条不可清除的记录。在本文中，为了简单起见，我们假设头中的所有记录都是可清除的。
请注意，我们知道日志尾部的每条记录都有一个唯一的键，因为重复的记录在上次清理中被删除了。只可能我们在 head 部分有一些记录，它们的键在日志中不是唯一的。为了更快地找到重复记录，Kafka 在 head 部分为记录创建了一个映射。回到我们的示例，偏移贴图结构如下所示:
如你所见，Kafka 创建了一个名为偏移贴图的结构，它为 head 部分中的每个关键点保存其对应的偏移。如果我们在头部有重复，卡夫卡使用最新的偏移。在上图中，键为 p6 的记录位于偏移量 5 处，p5 的最新偏移量为 7。现在，cleaner thread 检查日志中的每条记录，如果 offset map 中有任何记录具有相同的键，并且它的偏移量与 map 中的条目不同，就删除它(我们不想删除最新的记录)。
在压缩日志的清理过程中，不仅会删除重复的消息，而且 Kafka 还会删除值为 null 的记录。这些记录被称为墓碑。您可以通过设置delete.retention.ms 配置来延迟删除它们。通过设置此配置，Kafka 检查包含此记录的数据段的修改时间戳，如果修改时间早于配置值，则记录将被保留。
现在木头变得干净了。经过这个清洗过程，我们有了一个新的尾巴和一个新的头！为清理而扫描的最后一个偏移(在我们的例子中是旧头中的最后一个记录)是新尾的最后一个偏移。
Kafka 在数据目录的根目录下的一个名为cleaner-offset-checkpoint的文件中保存新头的起始偏移量。该文件用于日志的下一个清理周期。我们可以查看我们的主题检查点文件:
cat /var/lib/kafka/data/cleaner-offset-checkpoint
0
1
latest-product-price 0 6
如你所见，有三行。第一行是文件的版本(我认为是为了向后兼容)，第二行的值为 1，表示这一行之后有多少行(只有一行)，最后一行包含压缩日志主题的名称、分区号和这个分区的头偏移量。
结论
在本文中，我向您展示了什么是日志压缩主题，它们是如何存储的，以及 Kafka 如何定期清理它们。最后，我想指出的是，日志压缩非常适合于缓存场景，在这种场景中，您希望几乎实时地保存每条记录的最新值。假设您想在应用程序启动时构建缓存。您可以直接读取压缩的主题并构建您的缓存，因为 Kafka 按顺序读取消息，这比使用 SQL 数据库预热您的缓存要快得多。
你可以在 Martin Kleppmann 的文章《翻转数据库》中读到更多关于这种技术的内容。你可能还会发现我以前的文章使用 Kafka 和 Debezium 在 ASP.NET 内核中击败缓存失效很有用，它是这种技术的一个实现。
参考
https://github.com/apache/kafka/
https://thehoard . blog/how-kafkas-storage-internals-work-3a 29 b 02 e 026
https://medium . com/@ durgaswaroop/a-practical-introduction-to-Kafka-storage-internals-d5b 544 f 6925 f
https://kafka.apache.org/documentation/
数据线性化的对数变换基础无关紧要
Photo by Janus Clemmensen on Unsplash
解释为什么对数基数在线性化数据时没有显著影响的简单推导
此演示的代码可在此处找到:
今天一位同事问了我一个简单的问题:
“你如何找到最好的对数底数来线性变换你的数据？”
这实际上是一个棘手的问题，因为没有最好的对数基数来线性转换你的数据——不管对数的基数是多少，你取一个对数的事实都会使它线性化。
我的同事持怀疑态度，我想重温一下我的代数，所以让我们开始钻研数学吧！
前提
让我们假设你有指数数据。这意味着您的数据的形式类似于以下形式:
(1)
这意味着我们的数据是非线性的。线性数据可以说是我们可以建模的最佳数据形式，因为通过线性回归，我们可以通过查看目标变量的系数，直接量化每个特征对目标变量的影响。线性回归是最好的模型类型，它给人类一种直观和定量的感觉，即该模型认为我们的因变量是如何受我们的自变量影响的，例如，与深度神经网络的黑盒相比。
衍生物
因为我们知道这里的基数是 *e，*我们可以通过取两边的自然对数来线性化我们的数据(忽略常数 C₁):
(2)
现在，如果我们绘制 ln(y) vs. x ，我们会得到一条直线。这很简单，但是如果我们不知道我们力量的基础是什么呢？我们可以试着取两边的对数(基数为 10):
(3)
但是看起来还不是线性的。但是，如果我们引入对数幂法则呢？
(4)
但是 log(e) 是常数！因此，我们有:
(5)
这意味着我们以 10 为基数的对数仍然与 x 成正比，只是比例因子 C 不同，在本例中是原始基数 e 的对数。**
它看起来像什么？
我们也可以用一些 python 代码来可视化这一点！
import numpy as np
import matplotlib.pyplot as plt# Set up variables, x is 1-9 and y is e^x
x = list(np.linspace(1,10,100))
y= [np.exp(i) for i in x]# Plot the original variables - this is barebones plotting code, you # can find the more detailed plotting code on my github!plt.plot(x,y)# Plot log base 10 of y
plt.plot(x,np.log10(y))# Plot log base e of y
plt.plot(x,np.log(y))
They are both linear, even though the logarithms have different bases (base 10 vs base e)!
两个对数之间唯一变化的是 y 刻度，因为斜率略有不同！重要的是两者仍然与 x 成线性比例，因此在线性回归模型中具有相同的性能。
结论
综上:如果你有指数数据，可以做任意基数的对数变换，将数据线性化。如果你对领域知识的基础有直觉，那么就使用正确的基础——否则就没关系了。
附注:其他转换
如果您的数据是 x 的某个未知λ次方的稍微不同的形式，会怎么样？
(6)
在这种情况下， Box-Cox 变换将帮助您找到理想的幂来提升您的数据，以便将其线性化。我推荐用 Sci-py 的实现。
这一篇就讲到这里，感谢阅读！
逻辑理论—等价性
第四部分——使用真值表证明相等的陈述
Originally Published: https://www.setzeus.com/
本系列的第三部分以对真值表的初步理解而告终。我们看了两个独立的例子，映射出一个否定&一个暗示。然而，你可能已经意识到，虽然我们已经应用真值表来评估复合前提，我们还没有证明两个或更多陈述之间的任何等价性。
本文通过关注一个非常复杂的最终例子来扩展我们对真值表的了解。我们将直接比较两个复合前提*，以证明它们在逻辑上等价，而不是绘制出一个单个复合前提。具体来说，以下两种说法:*
陈述#1: (P ^ Q) - > P
语句#2: P - > (P ⌄ Q)
这一次我们将而不是为我们的两个前提提供进一步的上下文，就像我们在上一个例子中所做的那样——为什么？为了把不需要的观点讲清楚*。创建真值表的过程不*改变；我们严格关注源于每一个可能场景的个体真值&假值。
那么我们如何使用真值表来证明上述两者的等价性呢？
很简单，通过创建一个大规模的真值表来比较两个语句的最后两列。我们首先计算陈述#1 和陈述#2 各自的真值和假值；然后，之后，比较这些最终值，以证明(或证明)它们在逻辑上是等价的。
如果你想跳过，我们会在最后的真值表中清楚地看到这一点。不过现在，我们将从头开始。让我们来看看陈述 1。
构造语句#1
不要一次推导出(p ^ q)—>p 的所有场景，让我们把它分成两个主要部分。忽略语句的隐含部分(- > P) &我们只剩下一个更容易处理的语句。有了两个前提&和一个连接词(连词/and)，我们可以推断出我们的部分语句#1 表总共需要三列&四行:
Originally Published: https://www.setzeus.com/
如果你一直跟着做，上面的真值表是很简单的。我们在第一列和第二列简单地列出了两个前提(P，Q)所有可能的真值和假值，然后在第三列计算合取的真值和假值。回想一下，一个连词(P^Q)翻译成“P 和 Q 为真”正如你在上面看到的，当 P & Q 都为真时，这种情况只发生在第一行中。
无可否认，推导语句#1 的结果列，即蕴涵(P ^ Q) -> P，要复杂得多。下面，我做了一些非常规的事情，复制并重新定位了一个列(P ),以使计算的含义更加清晰:
Originally Published: https://www.setzeus.com/
过多地阅读专栏标题很容易让人迷失方向。尽管有多个变量，但这仍然是一个简单的暗示，就像我们在上一篇文章中回顾的那样。试着忽略列标题&从左到右，一次一行地分析隐含的意思。
第一行陈述了前因和后果都保持为真，这意味着整体蕴涵也确实为真。然而，在第二行到第四行中，antecedent 为 false。从我们之前的文章中，我们证实了在的任何蕴涵中，一个错误的前因总是导致一个正确的蕴涵；因此，第二到第四行中的隐含值也都是真。
我们走吧！有了所有真的最后一列，语句#1 现在被完全评估。假设的这一列为 true作为一个等式的左边，现在我们需要检查右边。证明语句#1 在逻辑上等价于语句#2 现在是一条清晰的前进道路。我们将计算语句#2 &的值，如果它们都为*真，*那么完美，我们有一个匹配；如果它们不同，那么我们可以明确地断定它们不等价。
构造语句#2
第二次就容易多了。我没有在多个步骤中写出表格，而是快进到完成的陈述#2 真值表:
Originally Published: https://www.setzeus.com/
从计算列数开始，比上一个示例中的语句#1 表少一列——不要让这使您迷惑，这次我们只是选择了而不是重写列§。此外，请注意我们从陈述#1 中重复的任何其他列？前提列§ & (Q)的两个都出现在语句# 1 &语句#2 中。最后但同样重要的是，看一下最后一列，这里的确实确认了所有真实值的另一列，与语句#1 中的完全相同。这最后一个观察回答了我们要回答的问题:由于两个语句在每个场景中都有相等的结果真/假值，我们可以无可争议地得出结论，它们在逻辑上是等价的。
最后
现在，这两种说法都推断出了它们最终的真值和假值，构建最终的真值表严格来说就是重新排列我们一路推导出来的列。如上所述，我们只需要写入我们的§和(Q)列*一次。*接下来，我们放入两个中间连接词，我们为两个语句分别计算了连接词(^) &析取词(⌄)。最后，我们有两个语句结果，用蓝色字母表示的两列，这是我们在比较语句逻辑上是否等价时唯一关心的两列。
通过第三个例子，我们最终应用了一个真值表来证明逻辑等价；当我们一点一点地建造它的时候，随着重复的曝光，这种需要会变得越来越少。理想情况下，当开始一个未来的问题时，我们至少会知道哪些列是必要的，以及如何安排它们——直接跳到上面抛光的真值表的结构。
如果你一直坚持到现在，那么恭喜你。虽然最后一个例子仅仅触及了逻辑证明的表面，但我们至少现在已经具备了在不久的将来进行更深入研究所需的基本理解。然而现在，在继续进行高级证明之前，我们将再次绕道一个名为 集合论 的迷人分支。照例评论，反馈&问题不胜感激。
来源
逻辑理论—真值表
第三部分——跨学科逻辑工具介绍
Originally Published On https://www.setzeus.com/
现在已经掌握了逻辑理论的原理以及基本符号，是时候探索逻辑中的等价的概念了。具体来说，是什么让两个复合前提相等？
两个复合前提 X & Y 是 逻辑等价的 如果，对于组成 X & Y 的原始前提的真值的每个赋值，陈述 X & Y 具有相同的真值。
这是一个难以接受的定义，但我们关心的是这个定义的应用。为了实现这一点，我们将遍历多个越来越复杂的例子。不过，首先，让我们绕道了解一下我们的 Excalibur 这是逻辑学家证明逻辑等价的最简单、最强大的工具之一: 真值表 。
真值表介绍
真值表是一种可视化工具，以具有行和列的图表的形式，显示复合前提的真或假。这是一种组织信息的方式，从提供的前提中列出所有可能的场景。让我们从最简单的例子开始，一个真值表描述了一个单一的前提操作:一个原始前提§的否定(~)
Originally Published On https://www.setzeus.com/
真值表总是从左向右读，第一列有一个原始前提。在上面的例子中，我们的原始前提( P )在第一列；而的结果前提( ~P 的 *)，*的后置否定，则构成了两列。
这里很容易想多了——不要忘记前提只是一个陈述，要么是真的，要么是假的。由于这个例子只有一个单个前提，我们只需要追踪两个结果；当 P 为真或为假时会产生两行。第一行描述，从左到右阅读，如果 P 为真，那么 P 的否定为假；第二行显示如果 P 已经为假，那么 P 的否定为真。
让我们继续看一个更复杂的例子，通过插入一个我们以前见过的连接词:蕴涵(->)。为了使这一点更容易理解，在建立真值表之前，让我们给我们的语句 P & Q 分配一些上下文:
灭霸打响指
问:50%的生物消失了
在看下面之前，考虑一下上面给出的细节。首先，因为我们有两个原始前提(P，Q)，我们知道我们至少需要两列；此外，我们应该为带有隐含连接词(P - > Q)的结果前提做准备，这将需要另一列。共有三列。
行呢？因为我们有两个前提，每个都可以为真或为假，为了说明所有可能的情况，我们总共需要 四个行*(P . S——从这个观察可以得出一个简洁的推论:说明 N 个前提的真值表需要 N 行)*。现在让我们把这张表画出来&确保它是可以理解的:
Originally Published On https://www.setzeus.com/
逐行回顾上面的真值表。第一行证实了灭霸打了他的手指§ & 50%的生物消失了(Q)。既然两个前提都成立，那么结果前提(蕴涵或条件)也成立:
Originally Published On https://www.setzeus.com/
第二排在理解上同样直接。这一次，P 仍然为真，然而 Q 现在为假。这里的解释是“灭霸掰着手指头，但是百分之五十的生物做到了 而不是 消失。”既然我们开始证明这种暗示的有效性，那么前面的陈述把整个前提完全错误是有道理的:
Originally Published On https://www.setzeus.com/
最后两行有点反直觉。这里有一个捷径:我们只需要看第一列，就可以证明这个蕴涵是真的。在第三行和第四行中，先行前提§为假— ，这是我们所需要知道的，不考虑前提 Q 的值，以确定蕴涵为真。
Originally Published On https://www.setzeus.com/
为什么一个错误的前因总是导致一个正确的暗示？因为在我们逻辑陈述的宇宙中，既然前因没有发生，就不可能排除所有可能导致 q 的可能场景。例如，第 3 行说“灭霸没有折断他的手指，然而 50%的生物消失了”。嗯，据我们所知，流星、自然灾害、外星人入侵或无数其他活动可能导致了那次灭绝——在这些场景的任何一个场景中，不管是哪一个，这种暗示都是真实的，因为我们仍然无法证明当他打响指时会发生什么。
证明等价性
真值表是光滑、方便的逻辑跟踪图，不仅出现在数学中，也出现在计算机科学、电子工程和哲学中。符号可能会因您从事的行业而异，但基本概念是相同的。它们是一种多用途的跨学科工具——然而我们仅仅触及了它们效用的表面。
现在有了真值表，是时候证明多个复合前提之间的等价性了。在本系列的下一篇文章中，我们将利用我们的复合知识来证明两个不同的复合前提，比如蕴涵式&反正，是相等的。
原载于
来源
作为非线性分类器的逻辑回归
逻辑回归传统上被用作线性分类器，即当类别可以在特征空间中由线性边界分开时。但是，如果我们碰巧对决策边界的形状有更好的想法，这是可以补救的……
逻辑回归是已知的，并被用作线性分类器。它用于在特征空间中提出一个超平面，以将属于某一类的观察结果与不属于该类的的所有其他观察结果分开。因此，判定边界是线性的。使用逻辑回归作为线性分类器的健壮和有效的实现是容易获得的(例如 scikit-learn)。
虽然逻辑回归对观察值做出了核心假设，如 IID(每个观察值都独立于其他观察值，并且它们都具有相同的概率分布)，但线性决策边界的使用不是之一。线性决策边界是出于简化的原因而使用的，遵循禅宗的原则——当有疑问时就简化。在我们怀疑决策边界是非线性的情况下，用非线性模型制定逻辑回归并评估我们能做得多好可能是有意义的。这就是这篇文章的内容。这是提纲。我们在这里浏览一些代码片段，但是可以从 github 下载重现结果的完整代码。
- 简要回顾似然函数的公式及其最大化。为了避免代数运算，我们在二维特征空间中坚持两个类。特征空间中的点[x，y]只能属于其中一类，函数 f (x，y；c) = 0 定义了决策边界，如下图 1 所示。
Figure 1. The decision boundary f(x,y; c) = 0 in feature space, separates the observations that belong to class A from those that do not belong to class A. c are the parameters to be estimated.
- 考虑一个决策边界，它可以用特征变量的多项式表示，但在权重/系数中是线性的。这种情况适合使用 scikit-learn 的 API 在(线性)框架内建模。
- 考虑一个不能用多项式表示的一般非线性决策边界。在这里，我们靠自己找到最大化似然函数的模型参数。scipy.optimize 模块中有很好的 API 可以帮助我们。
1.逻辑回归
逻辑回归是从一组连续和/或分类的观察值中预测(回归到)离散结果的一种练习。每个观察值都是独立的，并且一个观察值属于该类的概率是相同的！)描述该观察的特征的函数。考虑一组 n 的观测值[ x_i，y _ I；Z_i 其中 x_i，y_i 是第次次观测的特征值。如果第个个观测值属于该类，则 Z_i 等于 1，否则等于 0。获得 n 个这样的观测值的可能性仅仅是分别获得每个观测值的概率 p(x_i，y_i) 的乘积。
Equation 1. Define the likelihood & log-likelihood functions
比率 p/(1-p) (称为优势比率)沿着决策边界将是 1，即 p = 1-p = 0.5。如果我们定义一个函数 f(x，y；c) 为:
Equation 2. f is the functional form for log(odds ratio)
然后 f(x，y；c) =0 将是决策边界。 c 是函数 f 中的 m 参数。根据 f 的对数似然为:
Equation 3. Log-likelihood in terms of the data and parameters c. Ready to be maximized for the correct c
现在剩下要做的就是找到一组参数值 c ，使等式 3 中的 log(L) 最大化。我们可以应用最优化技术，或者为 c 求解耦合的非线性方程组md log(L)/DC = 0。
Equation 4 When LL is maximum, d log(L)/dc = 0 at that c value(s)
不管怎样，手里有了 c 就完成了 f 的定义，并允许我们找出特征空间中任意点的类别。这是坚果壳中的逻辑回归。
2.函数 f(x，y；c)
在等式 2 中， f 的任何函数都有效吗？肯定不会。随着 p 从 0 到 1， log(p/(1-p)) 从 -inf 到 +inf 。所以我们需要 f 在两个方向上都是无界的。只要我们确保 f 的范围从 -inf 到 +inf ，那么 f 的可能性是无穷无尽的。
2.1 f(x，y；c)
选择线性形式，例如
Equation 5. The default linear form used in general (e.g. scikit-learn)
肯定会起作用，这导致传统的逻辑回归可用于 scikit-learn，这也是逻辑回归被称为线性分类器的原因。
2.2 f(x，y；c)
一个简单的扩展方程 5 是使用更高次的多项式。一个二阶方程简单来说就是:
Equation 6. A simple extension that lends itself to embed nonlinearities for linear logistic regression
请注意，上述公式与线性情况相同， if 我们将 x、xy 和 y 视为问题的三个附加独立特征。这样做的动机是，我们可以直接应用 scikit-learn 的 API 来获得 c 的估计值。然而，我们在结果部分看到，以这种方式获得的 c 不如直接求解等式 4 获得的 c 。至于为什么，这有点神秘。但是在任何情况下，用 scipy 中的模块求解方程 4 都很容易。方程 4 所需的导数简单地表示为
2.3，y；c)
像 f(x，y；c) = sin(c_0x + c_2y) = 0 不起作用，因为它的范围有限。但是下面会。
Equation 7. A non-polynomial form for f
我们可以再次计算等式 4 的导数并求解其根，但有时更简单的方法是直接最大化等式 3 中的 log(L) ，这就是我们在接下来的仿真中要做的。
3.模拟
当 f(x，y；c) >小 _ 值)和那些不属于( Z = 1 当 f(x，y；c) < small_value)我们知道 f(x，y；c)。
训练和测试的数据各占 80/20。我们使用 scikit-learn 的 API 来比较我们想要比较的线性情况。
当选择求解 4 中的方程 c 或最大化方程 3 中的似然性时，我们使用 scipy.optimize 模块。下面代码片段中的 LL 和 dLLbydc 只是分别实现了等式 3 和 4。scipy 例程是为了最小化，所以我们在每种情况下否定符号，因为我们希望最大化可能性。
最后，我们求解 c ，从一个接近 0 的小的初始猜测开始。
4.结果
我们分别看多项式和其他一般情况。
4.1 f(x，y；c)= c0+C1 x+C2 y+C3 x+C4 x y+C5 y
我们得到一个椭圆作为以下 c 值的判定边界
我们基于上述函数为每个类生成 1000 个数据点，并对这些数据应用逻辑回归，以查看它对决策边界的预测。
pipenv run python ./polynomial.py 1000 2.25 -3.0 -18.0 35.0 -64.0 50.0
下面的图 2 显示了不同方法获得的数据分布和决策界限。绿色三角形描绘的轮廓完美地分隔了数据，即 F1 值为 1。它是通过求解 c 的等式 4 获得的。红线是通过传统的逻辑回归方法获得的，显然是一条直线。考虑到非线性，它尽力做到最好，F1 得分为 0.5……蓝色三角形描绘的轮廓很有趣。它是通过在扩充的特征空间上应用传统的逻辑回归而获得的。也就是说，“x”、“xy”和“y”被视为三个附加特征，从而线性化* f 用于 scikit-learn。它的 F1 得分为 0.89，还不错。但是不完全清楚为什么它应该比求解 dlog(LL)/dc = 0 更差。
Figure 2. The elliptical boundary separates the two classes. Solving dlog(LL)/dc = 0 yields the green contour that best matches the true decision boundary. The blue contour obtained solving the augmented linear form has over 10% error (why?) while the default application of logistic regression is 50% in error.
系数的实际值如下表 1 所示。c_0 换算成 1，这样我们就可以比较了。显然，求解 dlog(LL)/dc = 0 获得的系数最接近产生数据时使用的实际值。
Table 1. Coefficients obtained for the elliptic boundary
4.2 f(x，y；c) = c_0 + c_1 sin(c_2 x ) + c_3 y = 0
图 3 显示了运行模拟时的数据分布和预测边界
pipenv run python ./generic.py 1000 1.0 -1.0 1.0 -1.0
来自传统逻辑回归的直线获得 0.675 的 F1 分数，而 log(LL)的直接最小化获得满分。
Figure 3. The green line obtained by minimizing log(LL) perfectly separates the two classes.
下面的表 2 比较了模拟中获得的系数的值。有趣的是，C1 和 C2 的符号在实际值和最小化预测值之间是相反的。这很好，因为 sin(-k) = -sin(k)。
Table 2. Coefficients obtained for the generic case
5.结论
传统上使用逻辑回归来提出将特征空间分成类别的超平面。但是，如果我们怀疑决策边界是非线性的，我们可以通过尝试 logit 函数的一些非线性函数形式来获得更好的结果。求解模型参数可能更具挑战性，但 scipy 中的优化模块会有所帮助。
原载于 2019 年 3 月 13 日xplordat.com。
逻辑回归:自下而上的方法。
(特征工程思想——额外收获)
快速提问？
我们可以使用 RMSE 或 MSE 作为分类的评估指标吗？如果有，什么时候可以用？如果不是，为什么我们不能使用它？
我最后会透露答案。
对你有什么好处？
这篇文章将带你了解逻辑回归的细微差别，并让你熟悉特性工程的思想。我想澄清的第一件事是，逻辑回归的核心是回归，但它是一种分类算法。
我的告白:
这是我仅次于决策树的第二喜欢的算法，仅仅是因为它提供的简单性和健壮性。尽管有大量的分类算法，但这种过时的算法经受住了时间的考验，你可以看到它与神经网络的无缝集成是它非常特殊的地方。
我们来谈谈数据
为了达到目标，我挑选了一个数据集，它将允许我讨论特征工程。这是一个密码强度分类数据集。它有 3 个等级，0-弱，1-中等和 2-强。它只有两列，分别是密码和力量等级。该数据集有大约 700，000 条记录，因此可以肯定它不是一个玩具数据集。
免责声明！
这篇文章期望你能提供一些概率概念。请刷新 概率 (5 分钟)和 赔率 (8 分钟)。阅读 BetterExplained 的这些优秀文章: 指数函数直观指南& e 和 揭秘自然对数(ln) 。然后，复习一下这个 对指数函数和对数的 的简要总结。
OMG！我没有正确的数据。
这里可以关注代码。
陷阱 1:在密码数据集中可能有逗号(，)，并且您正在读取 CSV 文件，这是一个常见的分隔值文件。您可能会读取过多的逗号，并且可能会出现异常。所以，你可以通过设置 error_bad_lines = False 来避开这个问题。
Distribution of passwords based on their strength
我们拥有的数据是密码列表，熊猫将其识别为“对象类型”。但是，为了使用 ML 模型，我们需要数字特征。让我们从现有的文本中生成新的特征。这就是特征工程拯救我们的地方。特征工程通常基于直觉。做这件事没有固定的正确方法。通常，在行业中，人们倾向于依赖领域专业知识来识别新特性。根据我们的数据，你是领域专家。您知道理想密码的关键要素，我们将从这里开始。
创建的特征:
- 密码中的字符数
- 密码中的数字个数
- 密码中的字母数
- 密码中元音的数量
- 密码中辅音的数量(越多意味着密码的意义越小)
New features. We will have to see if they actually contribute something.
逻辑函数
因此，我们希望返回一个介于 0 和 1 之间的值，以确保我们实际上表示的是一个概率。为此，我们将利用逻辑功能。逻辑函数在数学上看起来像这样:
让我们来看看剧情
你可以看到为什么这是一个伟大的概率测量函数。y 值表示概率，范围仅在 0 和 1 之间。同样，对于 0 的 x 值，你得到 0.5 的概率，当你得到更多的正 x 值时，你得到更高的概率，更多的负 x 值时，你得到更低的概率。
利用我们的数据
好吧——这很好，但是我们到底该如何使用它呢？我们知道我们有五个属性——char _ count，numerics，alpha，元音，辅音——我们需要在逻辑函数中使用它们。我们可以做的一件显而易见的事情是:
其中 CC 是 char_count 的值，N 是 numerics 的值，A 是 alpha 的值，V 是元音的值，CO 是辅音的值。对于那些熟悉线性回归的人来说，这看起来很熟悉。基本上，我们假设 x 是我们的数据加上截距的线性组合。例如，假设我们有一个密码，其 char_count 为 8，numerics 为 4，alpha 为 4，元音为 1，辅音为 3。一些神谕告诉我们，𝛽0=1、𝛽1=2、𝛽2=4、𝛽3=6、𝛽4=8 和𝛽5=10.这意味着:
x = 1+(2 * 8)+(4 * 4)+(6 * 4)+(8 * 1)+(10 * 3)= 95
将此代入我们的物流功能，得出:
因此，我们认为具有这些特征的密码有 100%的可能性是中等(我取了第一行数据)
学问
好吧——有道理。但是是谁给了我们𝛽价值观呢？好问题！这就是机器学习中学习的用武之地:)。我们将了解我们的𝛽价值观。
步骤 1-定义你的成本函数
为了简单起见，我们将去掉一个类，使之成为一个二元分类问题。在本文的后面，我们将看到如何进行多重分类。现在，我把中等和弱混合作为一个类。如果你一直在研究机器学习，你可能会听到“成本函数”这个词。不过，在此之前，让我们思考一下。我们正在尝试选择𝛽值，以最大化正确分类密码的概率。这就是我们问题的定义。假设有人给了我们一些𝛽价值观，我们如何确定它们是否是好的价值观？我们在上面看到了如何获得一个例子的概率。现在想象一下，我们对所有的密码观察—所有的 700k 都这样做了。我们现在会有 70 万的概率分数。我们希望强密码的概率值接近 1，弱密码的概率值接近 0。
但是我们并不关心只得到一个观察值的正确概率，我们想要正确地分类我们所有的观察值。如果我们假设我们的数据是独立且同分布的，我们可以取所有我们单独计算的概率的乘积，这就是我们想要最大化的值。所以在数学中，如果我们定义逻辑函数和 x 为:
这可以简化为:
∏符号表示将产品归类为该密码进行观察。在这里，我们利用我们的数据被标记的事实，所以这被称为监督学习。另外，你会注意到，对于弱观测值，我们取 1 减去逻辑函数。这是因为我们试图找到一个最大化的值，由于弱观测值的概率应该接近于零，1 减去概率应该接近于 1。所以我们现在有了一个我们试图最大化的价值。通常，人们通过使其为负来将其转换为最小化:
注意:最小化负面和最大化正面是一样的。上面的公式将被称为我们的成本函数。
步骤 2 —渐变
现在我们有一个要最小化的值，但是我们实际上如何找到最小化我们的成本函数的𝛽β值呢？我们要不要试一堆？这似乎不是一个好主意…
这就是凸优化发挥作用的地方。我们知道物流成本函数是凸——请相信我。因为它是凸的，它有一个全局最小值，我们可以使用梯度下降收敛到这个最小值。
这是一个凸函数的图像:
source: utexas.edu
现在你可以想象，这条曲线是我们上面定义的成本函数，如果我们在曲线上选择一个点，然后沿着它下降到最小值，我们最终会达到最小值，这就是我们的目标。这里的是一个动画。这就是梯度下降背后的思想。
所以我们跟踪曲线的方法是，计算梯度，或者成本函数对每个𝛽.的一阶导数所以让我们做些数学计算。首先认识到，我们也可以将成本函数定义为:
这是因为当我们取对数时，我们的乘积变成了一个和。参见日志规则。如果我们定义𝑦𝑖在强观察时为 1，弱观察时为 0，那么我们只对强观察做 h(x ),对弱观察做 1 — h(x)。让我们对这个新版本的成本函数，对𝛽0.求导记住我们的𝛽0 在我们的 x 值里。所以记住 log(x)的导数是 1/𝑥，所以我们得到(对于每次观察):
使用商法则我们看到 h(x)的导数是:
x 对𝛽0 的导数是 1。综上所述，我们得到:
简化为:
引入负数和，我们得到对𝛽0 的偏导数是:
现在其他的偏导数就简单了。唯一的变化是𝑥𝑖的导数不再是 1。对𝛽1 来说是 CC𝑖，对𝛽2 来说是倪。所以𝛽1 的偏导数是:
为了𝛽2:
步骤 3 —梯度下降
现在我们有了梯度，我们可以使用梯度下降算法来找到最小化成本函数的𝛽s 值。梯度下降算法非常简单:
- 最初猜测您的𝛽值的任何值
- 重复直到收敛:
- 相对于𝛽𝑖的𝛽𝑖=𝛽𝑖−(𝛼∗梯度)
这里𝛼是我们的学习率。基本上，我们的成本曲线需要多大的步骤。我们现在做的是取当前的𝛽值，然后减去梯度的一部分。我们做减法是因为梯度是最大增加的方向，但是我们想要最大减少的方向，所以我们做减法。换句话说，我们在成本曲线上选取一个随机点，通过使用梯度的负值来检查我们需要向哪个方向靠近最小值，然后更新我们的𝛽值以向最小值靠近。重复直到收敛意味着不断更新我们的𝛽值，直到我们的成本值收敛或停止下降，这意味着我们已经达到了最小值。此外，同时更新所有𝛽值也很重要。这意味着您使用相同的先前的𝛽值来更新所有接下来的𝛽值。
梯度下降技巧
标准化变量:
- 这意味着每个变量减去平均值并除以标准差。
- 学习率:如果不收敛，学习率需要更小——但是收敛需要更长的时间。值得尝试的好值…，. 001，. 003，. 01，. 03，. 1，. 3，1，3，…
- 如果成本下降小于 10^−3，则声明收敛(这只是一个不错的建议)
- 绘制收敛图作为检查
遵循代码来实现。从零开始建立模型后，我们得到了 76%的准确率，我使用了 sklearn 软件包，得到了 99%的准确率，这相当不错。
高级优化
所以梯度下降是学习𝛽值的一种方法，但是还有其他的方法。基本上，这些是更高级的算法，我不会解释，但一旦你定义了你的成本函数和梯度，就可以很容易地在 Python 中运行。这些算法是:
我会把它留在这里让你去探索。我建议没有多少人知道这些先进的优化技术，而是他们很好地工作与梯度下降及其变种，因为他们做了相当不错的工作。
评估分类
现在，我们已经花了一些时间来理解逻辑回归背后的数学，让我们看看如何评估分类问题。我们现在可以拟合一个逻辑回归模型。逻辑回归也可以利用我们在线性回归文章中使用的正则化技术。它们的使用方式完全相同——我们通过向成本函数添加 L1 或 L2 范数来惩罚大系数。
请注意，我们不会像线性回归一样在这里演示超参数优化的 CV，但这仍然适用。事实上，我们在上一篇文章中讨论的几乎所有技术都可以在这里使用。我们将讨论新的话题。我们已经知道了前一篇文章中的训练和测试数据。让我们看看这个模型在测试数据上表现如何。
准确(性)
分类时最容易理解的度量之一是准确性:您正确预测的分数是多少。
Excellent!
准确性的一个缺点是，当你有不平衡的类时，这是一个相当可怕的度量。让我们来看看我们的阶级平衡:
Class Imbalance
因此，我们大约 88%的数据是负面的——因此，如果我们只是预测一切都是负面的，我们的准确率将是 88%左右！还不错。想象一下其他数据，其中负类只出现了 1%的时间——你可以通过预测所有事物都是正的来获得 99%的准确率。这种准确性可能感觉很好，但如果这是你的模型，它是非常没有价值的。
精确度、召回率和 F1
为了处理其中的一些问题，您经常会看到使用精度、召回率和 F1 来评估分类任务。
精确度是真阳性的数量除以所有的阳性预测(真阳性加上假阳性):基本上，当我们预测某事为阳性时，我们的预测有多精确。
Precision
回忆是真阳性的数量除以所有实际阳性的数量(真阳性加上假阴性):基本上，我们实际预测了所有阳性的多少。
Recall
这两个指标相互之间有所取舍。对于同一个模型，你可以以召回为代价提高精度，反之亦然。通常，您必须确定对于手头的问题哪个更重要。我们希望我们的预测是精确的还是高召回率的？
什么是真/假阳性/阴性？我们来看一个 混淆矩阵 。
Confusion matrix
这是我们测试预测的混淆矩阵。这些行就是真相。因此，第 0 行是实际的 0 标签，第 1 行是实际的 1 标签。列是预测。因此，单元格 0，0 计算我们正确得到的 0 类的数量。这被称为真正的否定(假设我们认为标签 0 是否定标签)。单元格 1，1 计算我们得到正确的 1 类的数量—真阳性。单元格 0，1 是我们的假阳性，单元格 1，0 是假阴性。
正如你可能知道的，混淆矩阵对于错误分析是非常有用的工具。在这个例子中，我们有非常对称的错误——两个类中都有四个错误。情况并不总是这样，通常会有两个以上的类，因此了解哪些类会让其他人感到困惑，对于改进模型非常有用。
有时候你只想在精确度和召回率之间取得平衡。如果这是你的目标， F1 是一个非常常见的指标。这是精确和回忆的调和平均值:
Harmonic mean of precision and recall
使用调和平均值是因为它强烈倾向于被平均的最小元素。如果 precision 或 recall 是 0，那么 F1 是 0。最好的 F1 成绩是 1。
精度/召回曲线
对于之前的所有指标，我们一直使用二元预测。但是真正的逻辑回归返回概率，这是非常有用的。
现在每个预测不是一个 1/0 的数字，我们有两个数字:被分类为 0 的概率和被分类为 1 的概率。这两个数字的和必须是 1，所以 1 的概率= 0 的概率。一件有用的事情是确保这些概率实际上似乎与类别相关联:
这张图向我们展示的是，对于我们的训练数据，我们的积极类确实有很高的概率是积极的，而我们的消极类有很高的概率是消极的。这基本上是我们设计算法的目的，所以这是可以预料的。检查验证集是否如此也是很好的。
Sklearn 的 predict 函数只是把概率最高的类作为预测类，这个方法还不错。然而，人们可能想要一个非常精确的模型来进行积极的预测。我们可以这样做的一个方法是，只有当概率超过 90%时，才称之为肯定的。为了理解不同截止值下精确度和召回率之间的权衡，我们可以绘制它们:
不错！当我们调整阈值时，我们可以清楚地看到 P 和 R 之间的权衡。事实上，我们可以找到最大化 F1 的截止点(注:我们在这里使用测试数据，这在实践中并不好。我们希望使用验证集来选择临界值)
受试者工作特征曲线
ROC 曲线是评估分类任务的另一种流行方法。它是一个图，y 轴是召回值(或真阳性率)，x 轴是假阳性率(假阳性除以所有实际阴性)。因此，它显示了您的模型如何在这两个值之间进行权衡。您的值越靠近图表的左上角越好。让我们来看看我们的测试数据:
ROC curve USP: It would work even with highly imbalanced data
显然，我们的模型做得非常好。您在右下角看到的 AUC(曲线下面积)分数是 ROC 曲线下的面积，1 是最好的值。我们有一个。如果您知道想要在 TPR 和 FPR 之间做出的权衡，该曲线也可用于选择阈值。
多类
逻辑回归也可以处理 2 个以上的类。有两种方法可以做到这一点:
- One-vs-rest 方法:使用这种方法，我们可以为每个类训练一个分类器，当该类是活动的时，正值，而对于所有其他类，负值。
- 交叉熵损失:我们也可以改变我们的损失函数来合并多个类。最常见的方法是使用交叉熵损失:
Cross entropy loss
这种损失采用正确类别的预测概率的负平均值。因此，优化它试图获得尽可能高的正确类别的预测概率。
解释逻辑回归
解释逻辑回归的系数比解释线性回归的系数要稍微复杂一些。让我们来看看我们的系数:
我们最正的系数是 char_counts 的值 5(粗略估计)。由于我们使用对数损失，我们需要通过 e^5 来转换这个值，即 147.52。也就是说，char_counts 每增加 1 个单位，密码是强密码的几率就会增加 14752%，因为我们现在已经转换为一个几率比。关于这个的更多细节，请看这个的帖子。
问题答案
我想用一个例子来说明。假设一个 2 类分类问题和 6 个实例。
假设，真实概率= [1，0，0，0，0]
案例 1: 预测概率= [0.2，0.16，0.16，0.16，0.16，0.16，0.16]
**情况二:**预测概率= [0.4，0.5，0.1，0，0，0]
情况 1 和情况 2 的均方误差分别为 0.128 和 0.1033 。
尽管如果 0.5 被认为是阈值，情况 1 在预测该实例的类时更正确，但是情况 1 中的损失高于情况 2 中的损失。这是为什么呢？MSE 是一种度量标准，用于测量预测值与实际值的标准偏差。而且它的应用可能仅限于二元分类。它仍然可能在特定的情况下使用。RMSE 是衡量分类器性能的一种方式。错误率(或错误分类的数量)是另一个因素。如果我们的目标是一个低错误率的分类器，RMSE 是不合适的，反之亦然。这解释了为什么我们不使用 MSE/RMSE 作为分类问题的评估指标。其他度量比 MSE 提供了更多的灵活性和对分类器性能的理解。
代号: Github
我的文章关于线性回归
如果你需要任何帮助，你可以联系: LinkedIn
参考资料:
Geron Aurelien 的动手机器学习
扫一扫
723
被折叠的 条评论
为什么被折叠?
到【灌水乐园】发言
