聊一聊登革热

March 22, 2020 zr9558 Leave a comment

近几个月，世界各地都爆发了新型冠状病毒（COVID-19），该病毒的传染性和致死率已经引起了许多国家的高度关注。其实，除了新型冠状病毒之外，在一些热带国家和地区，还有许多传染病值得大家注意防范，登革热就是人们需要关注的传染病之一。

登革热的英文名为 dengue fever，是一种由登革热病毒引起的传染病。其潜伏期的范围是 3 – 14 天（平均为 4 – 7 天），患者发病后会有发烧，头痛，关节痛，皮疹等诸多症状，一般情况下患者会在 2 – 7 天内痊愈。对于少数患者可能会出现病情加重的情况，甚至导致生命危险。登革热的种类分 1 – 4 型（DENV1-4），因此对于每个人来说，一生中最多获得四次登革热（但是个人觉得如果得了四次很可能也没法医治了）。对于年长者和体质虚弱的人，第一次获得登革热的时候将会比较严重，尤其需要注意防范登革热。

登革热作为热带地区的常见疾病之一，其传播途径却是相对有限的，比起新型冠状病毒可以人传人，登革热的传播途径其实只有蚊子。蚊子需要从登革热的感染者那里取血，然后经过 8 – 12 天的时间具备感染力，并且在其生命周期内都具有传染性。当一个健康的人被具有传染性的蚊子叮咬之后，就会在潜伏期之后（长达 3 – 14 天）发病，并且有一定比例的感染者不会有症状。经科学研究证明，登革热并不会直接人传人，因此与患者接触并不会传播疾病。传播登革热的蚊子主要是黑斑蚊（Aedes mosquitoes）,又称伊蚊，斑蚊，其中的白线斑蚊和埃及斑蚊就是登革熱的传播媒介。

要想阻断登革热的传播，从其传播媒介来看主要方法就是灭蚊。而灭蚊的方式可以使用蚊香，灭蚊液，防蚊的方法可以使用蚊帐等工具。同时每家每户都要保证家里的清洁卫生，例如时刻注意花瓶和水槽里面的积水，避免给蚊虫带来良好的生存繁殖环境。

moh_宣传材料_1

在新加坡的一些社区，随处可见 DENGUE 这样的宣传标语提醒人们注意预防登革热的传播。只有达到了足够的宣传力度，才能够时刻提醒人们的注意。

从 National Environment Agency（NEA）的报告数据中可以看出，在新加坡 2020 年 3 月份每天登革热的病例数大约是 28 – 64 人不等，每周大约有 300 – 400 人感染登革热。

相对于历年的登革热染病数据，2019 – 2020 年的登革热病例呈现增长的趋势，每周大约都有 300 – 400 人染病，相对于 2017 – 2018 年每周不到 100 人染病有了明显的上升。

新加坡面积不大，于是政府把相应社区的染病人数做了一个统计，让社区的人们注意安全，并且及时灭蚊。警告等级分三层：有十个以上案例就标记为红色，有十个以下案例就是黄色，没有案例就是绿色。

例如：Jurong West Street 91 的人数就较多，危险等级较高；Hougang Ave 8 的人数较少，危险等级较低。

nea_dengue_cases_4

不管是新型冠状病毒，还是热带地区常见的疾病，大家都要提高自我保护意识，做好防蚊等有效的措施。如果出现了身体不适，应该及早就医。

Singapore

新加坡首次出现两起新冠病毒死亡病例！

March 21, 2020 zr9558 Leave a comment

就在今天（2020 年 3 月 21 日），新加坡首次出现两起由于新冠病毒（COVID-19）引起的死亡病例。

从 MOH（Ministry of Health）的官网上可以看到，今天已经有两位感染新冠病毒（Coronavirus Disease 2019，简称为 COVID-19）的病人在新加坡离世。

第一名病人是一位 75 岁的女性新加坡公民，于 2020 年 2 月 23 日送入 NCID，并且当天确诊感染了 COVID-19。

第二名病人是一位 64 岁的男性印度尼西亚公民，于 2020 年 3 月 13 日送入 NCID，并且于 2020 年 3 月 14 日确诊感染了 COVID-19。

从新加坡政府的官网上可以看到：截止时间至 2020 年 3 月 21 日 12 点，现存确诊为 252 例（238 例稳定，14 例严重），死亡 2 例，治愈 131 例。

疾病爆发相应等级为橙色（orange）。

而新加坡政府也在其主页上呼吁大家注意卫生：

勤洗手，多用肥皂；
每天测量两次体温；
咳嗽或者打喷嚏的时候用纸巾遮住口鼻；
及时处理垃圾和食物，保持桌面整洁；
注意室内通风，保持厕所干净。

一旦发现自己生病，请带上口罩去医院，并且保证后续在家隔离。

而近三天，新加坡也有不少新增的病例。2020 年 3 月 20 日有 40 个病例，其中 30 位有欧洲，北美，亚洲的旅游史，7 位联系上了以前的病例，3 位没有找到传染病联系。

除此之外，从丁香医生的小程序可以看出，亚洲很多国家的疫情其实不容乐观。

最后，无论在国内还是国外，大家都要注意安全。抗击疫情，人人有责。

口号_1

智能运维

FluxRank: 如何快速地进行机器故障定位

March 21, 2020 zr9558 Leave a comment

在运维领域，服务侧的异常会由多方面的原因造成，有的时候是因为网络的抖动，有的时候是因为机器的故障，有的时候甚至是因为人为的变更。本篇博客会介绍一种机器异常定位的方法，论文是来自于清华 Netman 实验室的《FluxRank：A Widely-Deployable Framework to Automatically Localizting Root Cause Machines for Software Service Failure Mitigation》。本篇论文主要介绍了如何从服务的故障定位到局部异常的机器，也就是说在发现服务故障的同时，进一步推断出是由哪些机器出现问题而导致的。

通常来说，在服务异常（例如服务的耗时长，失败数上涨）的时候，需要运维人员通过历史上的经验迅速定位到是哪个业务，哪个模块，甚至哪台服务器出现了故障。而人工定位的速度总是会出现瓶颈的，无论对模块的判断，还是机器的判断，都依赖于人工所积累的经验。而每个人的经验却各不相同，并且经验的传承也需要一定的时间成本。那么如何基于人工运维的经验来构建模型，进一步地提升异常定位的速度就是智能运维的关键之处之一。

对于一条业务指标（时间序列）而言，大多数情况下是处于正常的状态（normal）。但是如果出现了错误的变更，发布了错误的程序，或者服务器突然出现了故障，都会导致业务指标出现变化，就从正常（normal）变成异常（abnormal）。这个时候就会出现一个故障的开始时间，也就是 failure start time $T_{f}$ ，这个时间戳是运维领域非常重要的时间戳，它由异常检测（anomaly detection）产生，无论在告警收敛（alarm convergence）还是根因分析（root cause analysis）都非常依赖这个时间戳。而另外一个时间戳虽然没有故障开始时间那么重要，但是也有着其实用价值，那就是缓和开始时间（mitigation start time），它表示故障虽然还没有恢复，但是出于稍微平稳的走势，并没有持续恶化。在出现了故障之后，通常都会发送相应的告警给运维人员，那么在发送告警的时候，如果将异常定位的结果随之带出，则会大大减少运维人员排障的时间。在故障缓和的时间内，运维人员通常需要进行必要的操作来排查故障，例如切换流量（switch Traffic），回滚版本（Rollback Version），重启实例（Restart Instances），下线机器等操作。除此之外，为了定位问题（Root Cause Analysis），运维人员需要分析源码（Code Analysis），查看日志（Log Analysis）等一系列操作。如果能够将这一系列操作融入相应的机器学习模块中，将会节省运维人员大量的排障时间。

贝叶斯网络

通常来说，故障定位也称为根因分析或者根源分析（Root Cause Analysis），都是为了排查产生这次故障的原因。在机器学习领域，为了进行因果分析（Causal Analysis），则需要使用相应的模型来进行建模。其中较为经典的统计分析方法则是贝叶斯分析法，其中的贝叶斯网络（Bayesian Network）则是经典模型之一。下面来看一个简单的例子。

假设降雨（Rain）的概率是 0.2，不降雨的概率是 0.8；而洒水器（Sprinkler）是否开启会受到降雨的影响，其条件概率与下图所示。而降雨或者洒水器都会导致草湿润（Grass Wet），其概率分布如下图所示。那么可以问如下问题：

如果草已经湿润，求降雨的概率是多少？
如果草已经湿润，求没有降雨且洒水器开启的概率是多少？

而这一类的问题可以通过贝叶斯公式来进行解答。从表格来看：

从 Rain 的表格可得： $P(R=T)=0.2, P(R=F)=0.8$ 。

从 Rain 和 Sprinkler 的表格可得： $P(S=T|R=F)=0.4, P(S=F|R=F)=0.6$ ， $P(S=T|R=T)=0.01, P(S=F|R=T)=0.99$ 。

针对问题 1，需要计算条件概率 $P(R=T|W=T)$ 。从 Bayes 公式可以得到： $P(R=T|W=T) = P(R=T,W=T)/P(W=T)$ 。分别计算分子分母即可：

$P(R=T,W=T)=P(R=T,S=T,W=T)+P(R=T,S=F,W=T)$

$= P(W=T|R=T,S=T)P(S=T|R=T)P(R=T) + P(W=T|R=T,S=F)P(S=F|R=T)P(R=T)$

$= 0.99*0.01*0.2+0.8*0.99*0.2=0.16038$ ，

$P(W=T)=P(W=T,S=T,R=T)+P(W=T,S=F,R=T)+P(W=T,S=T,R=F)+P(W=T,S=F,R=F)$

$= 0.99*0.01*0.2+0.8*0.99*0.2+0.9*0.4*0.8+0.0*0.6*0.8=0.44838,$

那么如果草已经湿润，求降雨的概率是 $P(R=T|W=T)=P(R=T,W=T)/P(W=T)=0.16038/0.44838=0.3577.$

另外一个题目可以用类似的方法进行求解，在此不再赘述。

虽然贝叶斯算法能够计算出条件概率，例如本次故障是由哪些原因导致的，但是这个需要长期收集数据，需要对历史数据进行积累，才能通过人工或者统计的方法得到以上表格的条件概率。但是在实际的环境中是较难获取这些数据的，需要大数据平台的支持，因此需要探索其他的解决方案。

FluxRank

在本论文中，为了克服贝叶斯网络模型中的一些问题，针对子机异常定位的场景，设计了一套技术方案，作者们称之为 FluxRank。

FluxRank 这一模块的触发需要服务指标（Service KPI）的异常，因此需要对服务指标（Service KPI）进行异常检测。这里的服务指标通常指的是业务指标，包括某块 APP 的在线人数，某个接口的成功率，某个视频网站的卡顿数等指标。当服务指标出现了异常的时候，就启动 FluxRank 模块进行异常机器定位。

如果按照人工处理的流程来看，分成几个步骤：

异常检测部分：通过设定阈值或者某个简单的规则来进行异常检测，包括服务的 KPI（Service KPI）和机器的 KPI（machine KPIs）；
手工检查异常的时间段，并且查看在异常的时间段内发生了什么情况；
运维人员根据自身的业务经验来对机器的故障程度做人工排序；
运维人员根据自身的业务经验来对故障进行处理，并且人工给出处理方案。

那么 FluxRank 所面临的挑战就有以下几点：

如何衡量海量 KPIs 的变化程度？在这里不仅有服务的 KPIs，还有机器的 KPIs。而机器的 KPIs 包括内存，硬盘，IO，CPU等诸多固定的指标，那么如何对这些海量的 KPI 曲线进行变化程度的衡量，为后续的指标排序做准备就成为了一个难点；
如何对 KPIs 进行异常性或者重要性的聚类，让运维人员能够一眼看出每个聚簇的差异或者异常程度？
如何对 KPIs 聚类的结果进行排序？

为了解决以上的问题，FluxRank 的框架有以下几个贡献点：

基于 Kenel Density Estimation 用于衡量海量 KPIs 在某一个时间段的变化程度和异常程度；
基于上一步生成的异常程度，对诸多机器所形成的特征使用距离公式或者相似度公式，然后使用 DBSCAN 聚类算法来对机器进行聚类；
在排序部分，对上一步的机器聚类结果进行排序；

Change Quantification

首先，来看一下 Change Quantification 是怎么样做出来的。这里的 Change Quantification 使用与衡量机器 KPIs 的变化程度，称之为 change degree。Change degree 可以用于 CPU，内存，IO 等诸多机器指标。为了达到衡量变化程度，需要一个非常重要的信息，那就是变化的开始时间，change start time，也就是说在哪个时刻时间序列开始出现了变化。于是在 Change Quantification 部分，就分成两部分：（1）用 absolute derivative 或者 CUSUM 算法获得变化开始时间（change start time）；（2）用 Kernel Density Estimation（KDE）来计算变化程度（change degree）。

正如上图所示，针对服务 KPIs（ervice KPIs），存在两个关键的时间点，那就是失败开始时间（Failure Start Time） $T_{f}$ 和缓和开始时间（Mitigation Start Time） $T_{m}$ 。在失败开始时间 $T_{f}$ 之前，可能有的机器已经出现了故障，因此变化开始时间（Change Start Time） $T_{c}$ 小于或者等于 $T_{f}$ 。通常情况下，一个或者多个机器故障会在半小时（30 mins）甚至更短的时间内引发服务故障，因此，只需要假设 $w_{1}=30$ 即可。关键时间点的排序为 $T_{f}-w_{1}<T_{c}\leq T_{f}<T_{m}$ 。

对于服务 KPIs 的异常检测，FluxRank 中提到了两种方法：分别是 absolute derivative 和 CUSUM 方法。

absolute derivative 方法：个人理解就是对时间序列进行一阶差分操作，然后对一阶差分来做时间序列异常检测，例如 3-sigma 等方法，一旦有明显的变化，就说明当前的时间点出现了突增或者突降；与该方法比较类似的一种方法是：MAD（Median Absolute Deviation）。对于一条时间序列 $X=[x_{1},\cdots,x_{n}]$ 而言，MAD 定义为 $MAD = median_{1\leq i\leq n}(|x_{i}-median(X)|)$ ，而每个点的异常程度可以定义为： $s_{i}=(x_{i}-median(X))/MAD = (x_{i}-median(X))/median_{1\leq i\leq n}(|x_{i}-median(X)|).$ 当 $s_{i}$ 较大或者较小的时候，表示上涨或者下降的异常程度。通过设置相应的阈值，同样可以获得时间序列的异常开始时间。
CUSUM 算法也是用于时间序列异常检测的。对于一条时间序列 $X=[x_{1},x_{2},\cdots,x_{n}]$ ，可以预估它的目标值（target value） $\mu_{0}$ ，通常可以用均值来估计，也需要计算出这条时间序列的标准差 $\sigma$ 。通常设定 $\mu_{1}=\mu_{0}+\delta\sigma$ ， $K=\delta\sigma/2=|\mu_{1}-\mu_{0}|/2$ 。而 Tabular CUSUM 指的是迭代公式 $C_{i}^{+}=\max[0,x_{i}-(\mu_{0}+K)+C_{i-1}^{+}]$ ， $C_{i}^{-}=\max[0,(\mu_{0}-K)-x_{i}+C_{i-1}^{-}]$ ，初始值是 $C_{0}^{+}=C_{0}^{-}=0$ 。当累计偏差 $C_{i}^{+}$ 或者 $C_{i}^{-}$ 大于 $H=5\sigma$ 的时候，表示 $x_{i}$ 出现了异常，也就是 out of control。通过这个值，可以获得时间序列开始异常的时间。

从论文的描述来看，作者是使用 absolute derivative 来做异常检测的，并且定位其异常开始时间的准确率较高。

Change Degree

其次，我们来看一下变化程度（Change Degree）是怎么计算出来的，通过之前的计算，我们已经可以获得一些关键的时间戳，例如 $T_{f}, T_{c}, T_{m}$ 等时间戳。根据变化开始时间（change start time） $T_{c}$ ，同样需要设置一个窗口值 $w_{2}$ ，例如 60 分钟（1 小时）。可以从两个时间段获取数据，正常时间段 $[T_{c}-w_{2},T_{c})$ ，异常时间段 $[T_{c},T_{m}]$ ，分别获取到数据 $\{x_{i}\}$ 和 $\{x_{j}\}$ ，前者是在变化开始时间之前的数据点，后者是在变化开始之后的数据点。于是，作者们通过概率值来计算变化程度 $P(\{x_{j}\}|\{x_{i}\})$ ，意思就是计算一个条件概率，在观察到 $\{x_{i}\}$ 之后，得到 $\{x_{j}\}$ 的概率值。

为了计算以上概率值，需要简化模型，因此这里需要假设 $\{x_{j}\}$ 是独立同分布（iid）的，于是 $P(\{x_{j}\}|\{x_{i}\})=\prod_{j=1}^{\ell}P(x_{j}|\{x_{i}\})$ ，在这里 $\ell$ 表示集合 $\{x_{j}\}$ 的元素个数。为了分别得到其上涨和下降到概率，则需要计算：

$P_{o}(\{x_{j}\}|\{x_{i}\}) = \prod_{j=1}^{\ell}P(X\geq x_{j}|\{x_{i}\})$ ,

$P_{u}(\{x_{j}\}|\{x_{i}\}) = \prod_{j=1}^{\ell}P(X\leq x_{j}|\{x_{i}\})$ ,

其中 $P_{o}(\{x_{j}\}|\{x_{i}\})$ 表示上涨的程度， $P_{u}(\{x_{j}\}|\{x_{i}\})$ 表示下降的程度。如果不想处理连乘的话，则需要处理连加：

$o=-\frac{1}{\ell}\sum_{j=1}^{\ell}\ln P(X\geq x_{j}|\{x_{i}\})$ ,

$u =-\frac{1}{\ell}\sum_{j=1}^{\ell}\ln P(X\leq x_{j}|\{x_{i}\})$ .

在这里，作者们使用了三种概率分布函数，分别是 Beta 分布（Beta distribution），泊松分布（Poisson distribution），高斯分布（Gaussian distribution）。

Beta 分布的概率密度函数（probabilisty density function）是 $f(x;\alpha,\beta) = x^{\alpha-1}(1-x)^{\beta-1}/B(\alpha,\beta)$ ，其中 $B(\alpha,\beta)=\Gamma(\alpha)\Gamma(\beta)/\Gamma(\alpha+\beta)$ 。在机器 KPIs 中，CPU 等指标可以用 Beta 分布；

泊松分布的概率密度函数是 $f(x;\lambda)=\lambda^{x}e^{-\lambda}/x!$ ，在机器 KPIs 中，SYS_OOM 用于衡量超出内存的频率，可以用泊松分布来做。

高斯分布的概率密度函数 $f(x;\mu,\sigma) = e^{-(x-\mu)^{2}/2\sigma^{2}}/(\sqrt{2\pi}\sigma)$ 。

根据论文中的陈述，机器 KPIs 分别适用于以下概率分布：

通过以上公式，可以计算出每一个机器的每一个指标的 $o$ 和 $u$ 两个值。

Digest Distillation

再来看一下 Digest Distillation 部分，在此部分需要对机器的 KPIs 进行聚类操作；那么就需要构造特征向量和距离函数，再加上聚类算法即可获得结果。

每一个机器的特征向量是由之前计算的 Change Degree 形成的，由于每台机器的 KPIs 都是一样的，因此可以对它们的 KPIs 的 change degree 进行排列。假设每台机器有 $k$ 个 KPIs，那么这台机器所对应的向量就是 $(o_{0},u_{0},\cdots,o_{k},u_{k})$ 。

在描述向量的相似性方面，可以使用相关性的系数，包括 Pearson 系数，Kendall tau 系数，Spearman 系数。对于两条时间序列而言， $X=[x_{1},\cdots,x_{n}]$ 和 $Y=[y_{1},\cdots,y_{n}]$ ，

Pearson 系数指的是： $\rho_{X,Y}=\sum_{i=1}^{n}(x_{i}-\overline{x})\cdot(y_{i}-\overline{y})/\sqrt{\sum_{i=1}^{n}(x_{i}-\overline{x})^{2}\cdot\sum_{i=1}^{n}(y_{i}-\overline{y})^{2}},$ 其中 $\overline{x}=\sum_{i=1}^{n}x_{i}/n$ ， $\overline{y}=\sum_{i=1}^{n}y_{i}/n$ 。

Kendall tau 系数指的是：如果 ( $x_{i}>x_{j}$ 且 $y_{i}>y_{j}$ ) 或者 ( $x_{i}<x_{j}$ 且 $y_{i}<y_{j}$ )，那么称之为 concordant；如果 ( $x_{i}<x_{j}$ 且 $y_{i}>y_{j}$ ) 或者 ( $x_{i}>x_{j}$ 且 $y_{i}<y_{j}$ )，称之为 discordant；如果 $x_{i}=x_{j}$ 或者 $y_{i}=y_{j}$ ，则既不是 concordant，也不是 discordant。那么 Kendall tau 定义为 $[\text{(number of concordant pairs)}-\text{(number of disordant paris)}] / [n(n-1)/2]$

Spearman 系数指的是：通过原始序列变成秩次变量（rank）（从大到小降序排列即可）， $x_{i}$ 将会对应到 $x_{i}'$ ，后者表示 $x_{i}$ 在从大到小排序之后的序列 $\{x_{i}\}_{1\leq i\leq n}$ 的位置，称之为秩次（rank），得到序列 $X'=[x_{1}',\cdots,x_{n}']$ 。对原始序列 $Y=[y_{1},\cdots,y_{n}]$ 作同样的操作，得到 $Y'=[y_{1}',\cdots,y_{n}']$ 。一个相同的值在一列数据中必须有相同的秩次，那么在计算中采用的秩次就是数值在按从大到小排列时所在位置的平均值。如果没有相同的 rank，那么使用公式 $r_{s} = 1-6\sum_{i=1}^{n}d_{i}^{2}/(n(n^{2}-1))$ 进行计算，其中 $d_{i}=x_{i}'-y_{i}'$ ；如果存在相同的秩次，则对 $X'=[x_{1}',\cdots,x_{n}']$ 和 $Y'=[y_{1}',\cdots,y_{n}']$ 来做 Pearson 系数即可，也就是 $\rho_{X',Y'}$ 。

通过作者们的实验，说明 Pearson 系数在这个数据集上效果最佳。在聚类算法的场景下，作者们同样对比了 KMeans，Gaussian Mixture，Hierarchical Clustering，DBSCAN 算法的效果，最后使用了 DBSCAN 的聚类算法。每一个聚类的结果，作者称之为一个 digest，也就是下图的 M1，M2 等聚类结果。

Digest Ranking

最后，就是对聚类结果的排序工作。通过观察会发现：

变化开始时间（change start time） $T_{c}$ 会在失败发生时间 $T_{f}$ 之前；
不同的故障机器 KPIs 的 change start time 是非常接近的；
故障机器的一些 KPIs 的 change degree 是非常大的；
故障机器的占比是与故障原因相关的，故障机器越多说明故障越大；

在同一个模块下，如果出现故障机器的占比较大，那么故障将集中于这个模块下，可以通过 ratio 这个指标进行排序工作。

实验数据

在 FluxRank 论文中，作者们收集了 70 个真实的案例，然后根据实验效果获得了结果。

在标记的时候，除了标记异常机器（Root Cause Machines，简称为 RCM）之外，也需要标记相关的指标（Relevant KPI，简称为 RK）。Root Cause Digest（简称为 RCD）把包括两个部分，不仅包括 RCM 的一个聚类结果，还包括聚类结果中的 top-five KPIs。

通过对 FluxRank 进行实验，可以得到如下实验数据：

其中 Recall@K 指的是： $Recall@K=\text{\# of cases whose top-k digests contain RCDs}/ \text{\# of all cases},$ 或者 $Recall@K=\text{\# of cases whose top-k machines contain RCMs}/\text{\# of all cases}.$

参考资料

FluxRank: A Widely-Deployable Framework to Automatically Localizing Root Cause Machines for Software Service Failure Mitigation，Ping Liu，Yu Chen，Xiaohui Nie，Jing Zhu，Shenglin Zhang，Kaixin Sui，Ming Zhang，Dan Pei，ISSRE 2019， Berlin, Germany, Oct 28-31, 2019。
Introduction to Statistical Quality Control，6th edition，Douglas C.Montgomery。
Bayesian Network：https://en.wikipedia.org/wiki/Bayesian_network

PHD的生涯

在新加坡的这五年—生活篇（五）

March 14, 2020 zr9558

在新加坡，民众运动健身的氛围极其浓厚。无论在周末的公园里，在工作日傍晚的的林荫跑道上，还是在各个区域都拥有的游泳池里，随处可见运动的人群。而在校园里，很多学生都会选择在操场上踢球，玩飞盘，搞活动。NUS 的足球场（Sports Field）就在数学系的旁边（S17），对于学生们而言，下楼的时候除了可以在售货机购买罐装饮料之外，还可以在操场旁边放空自己的大脑。当时的售货机在 NUS 是随处可见的，不过在移动支付尚未流行的年代，只能用硬币去购买各种饮料。整体来看，在校园里生活的话大家还是用现金居多，毕竟跟外出购物的大额消费不一样，日常买饭和坐公交车出行其实都可以用现金解决问题。结果硬币越来越多，钱包越来越重，但是钱却越来越少。于是笔者为了把零散的硬币（一角，两角，五角）这种换成一元硬币，除了去门店找人兑换之外，随处可见的售货机就是一个还不错的选择。

通过 Google 地图来搜索 Swimming Pool，会发现游泳池遍布新加坡全岛。除了最著名的滨海湾金沙酒店的无边际游泳池外，政府还在全岛修建了无数的游泳池供市民健身娱乐。当把地图聚焦于新加坡国立大学，将会看到其实 NUS 就有好几处游泳池，除了 UTown 的屋顶游泳池与金沙酒店的极其相似之外，在大学运动中心（University Sports Centre）也有一个标准的游泳池。而且这个游泳池是在 2015 年才开始重建，2018 年左右才开始对外开放。相比于原先的游泳池，其基础设施已经大大的改善了。

作为一个博士生，除了从办公室步行去学校的游泳池之外，更多的时候则是在周末的时候从宿舍坐车去游泳池。而 NUS 其实有着自己的穿梭巴士（Shuttle Bus），同时每一辆穿梭巴士其实都有着自己的路线。而周末从 PGP 坐车去游泳池的时候，可以选择早上 9:00 出发，坐上 A1 的即可在十分钟内到达游泳池。不过周末的巴士半小时才会行驶一趟，一旦错过了，只能够等待半小时以上了，否则就只能够步行前往目的地。

从 NUS 的巴士线路图上可以看出，NUS 不仅无缝对接着地铁站 Kent Ridge MRT，也对接着 Kent Ridge Bus Terminal。因此无论从学校的任何地方，都可以通过 Shuttle Bus 到达地铁站或者公交站。不过话说回来，笔者刚到 NUS 的时候是 2010 年，当时还没有开通 Kent Ridge MRT，要坐地铁只能够去 Clementi 或者 Buona Vista，而且只能坐 95 路或者 96 路公交车出行。到了 2011 年的时候，随着修建了 6 年的 Kent Ridge MRT Station 建成，NUS 的学生们就可以很方便的出行了。

NUS 的 Shuttle Bus 是免费给所有人乘坐的，包括老师，学生，外来的工作人员等。除了免费这一个福利之外，NUS 的 Shuttle Bus 也别具特色，学校会将 LOGO 印在 Bus 上，可以让人一眼认出这个是 NUS 的 Bus。而 Bus 的大小是按照线路来分的，通常来说 A1，A2 的车型相对偏大，而 B1，B2，C 车的车型相对偏小，这也是根据是否经过宿舍区，和路线上人数的多少而定的。

NUSShuttleBus_3

在 NUS 里面，除了学校的穿梭巴士之外，还有很多其余的穿梭巴士，其中比较有名的就是 West Coast Plaza 的 Shuttle Bus 了。在 UTown 和 Kent Ridge MRT Station 还没有建好的时候，West Coast Plaza 就是 NUS 学生经常去购物的地点了。因此，免费的穿梭巴士就成为了首选。但是 West Coast Plaza 的穿梭巴士则小了许多，跟香港的小巴能够装载的人数差不多。除此之外，在 PGP 附近，笔者好像也见过 Science Park Shuttle，只不过从来就没有坐过，也并不知道这辆车的目的地是哪里。

NUSShuttleBus_5 — West Coast Plaza Shuttle Bus

在新加坡，除了随处可见的游泳池之外，新加坡最著名的游泳池应该来自于金沙酒店了。金沙酒店位于新加坡的滨海湾，属于新加坡的核心地带，也是旅游人士必去的打卡地点。金沙酒店除了提供住宿之外，还提供饮食，娱乐等诸多活动。每到夜晚，滨海湾的夜景总是让人觉得美不胜收。

屋顶上的游泳池则是金沙酒店的亮点之一，有机会的话自然会上去游玩一番。在 2013 年 11 月 16 日，笔者和同学 XU zhe 同去金沙酒店游泳池游玩。游泳池位于金沙酒店 55 楼层高的塔楼楼顶，长度大约 198 米。当时为了同时欣赏白天和晚上的景色，选择了从下午 4 点待到晚上 8 点。对于想来新加坡游玩的旅行者，如果想要入住金沙酒店的话，其实可以提前 45 天左右预定，将会获得不少的折扣。只要入住了金沙酒店，就可以去上面的游泳池进行观光和游泳。

对于宾客而言，在进游泳池的时候，凭借房卡里面的游泳券就可以入场，工作人员会在你的手上缠着一个防水胶带，表示你今天随便出入游泳池。进场之后，可以向工作人员要一块免费的浴巾，然后就去卫生间里面换衣服。换了衣服之后，就可以出来在外面找一个免费的躺椅，把自己的随身物品搁置在上面，便可以在无边游泳池游玩了。为了保证大家能够观光和游玩，游泳池的水深才 1.2 米，成年人可以站在游泳池里走到边缘俯瞰整个滨海湾。

虽然金沙酒店的游泳池的边缘看上去是无边的，其实在游泳池的外侧有栏杆和水槽，不仅可以接住游泳池的水，还能够进行一定的安全防护。滨海湾的一栋栋高楼就是各大银行的办公地点，还有国家大剧院，体育场等诸多设施。泳池的两侧也安装了照明设施，每当夜幕降临，游泳池的灯光将会打开，不过此时游泳观光的人群就相对少了许多。

NUS Utown 的游泳池虽然没有金沙酒店的游泳池那么高端，但是也是 NUS 的一大亮点之一。UTown 是 NUS 于 2011 年修建的大学城，承担了教学，住宿，学习，娱乐，健身于一身的诸多功能。在 UTown 建成之前，Prince George Park 其实只提供给本科生住宿，而 master 或者 PHD 想住学校宿舍的话，就只能够选择 BoonLay 的宿舍或者 Common Wealth 两处的宿舍。因为宿舍实在是不够用，必须要优先解决本科生的住宿问题。当 UTown 建成了之后，有不少的博士生都选择搬进去，而搬进去的时候，可以选择自己的室友。如果是熟悉的朋友自然是好，如果是没有熟悉的朋友就只能够跟其他陌生人一起居住了，不过就算是陌生人也是 NUS 的学生，因此整体来看也不会出现太大的问题。

UTown 的游泳池就位于 Stephen Ready Centre 的第三层，它对学生免费的。因此只要是 NUS 的教职工，就可以通过证件在开放时间内进入游泳池。一般到了午后或者周末，游泳池里面总是有很多学生在游泳或者晒太阳。相对于 Swimming Pool，这里更像是一个休闲的地方。无论是学生还是老师，觉得科研或者学习疲惫了，就可以来游泳池运动并且放松一下。

随着夜色的来临，其实 UTown 的游泳池也有着自身独特的魅力。不仅可以看到整个 UTown 的夜景，还能够在放松之余思考人生。除了科研，吃喝玩乐之外，思考人生则是博士生日常生活的重要组成部分。

在 UTown 里面，除了游泳池之类的健身场所之外，还有一个 7*24 小时开设的星巴克。该门店位于 Education Resource Centre（ERC），能够容纳 300 多个人。每个学期到了 Recess Week 或者 Reading Week 的时候，星巴克就是人山人海。不少学生都会选择在学校努力学习，奋笔疾书以应对即将来临的考试。而此刻星巴克外部的桌位通常来说都会被学生们坐满，学生们会在外部的长凳子上进行小组学习（Group Study），通过小组讨论或者自学的方式来复习本学期相应的知识点。

其实 NUS 的食堂还是相对较多的，在 Engineering，Science，Business 等院系都有自己的食堂。无论是学校内的人，还是外部的工作人士，都可以选择到这些食堂吃饭。在大学城 UTown 的 Town Plaza 也有食堂，无论是中国美食还是其他国家的食物，基本上在这里都能够找到。对于数学系的博士生而言，每天到了午餐或者晚餐的时候，除了在 Science 吃饭之外，还可以选择到 UTown 或者 Prince George Park 吃饭。而每到周末或者节假日，也就只有这两个食堂开门了，其他食堂（例如 Science Canteen 等食堂）在节假日期间是不会对外开放的。在博士生日常思考人生的同时，思考每一顿饭吃什么也是一件非常重要的事情。

UTown 的 Graduate Residence 则是为学生们开设的宿舍区，每个屋子里面有四个小房间，每个房间里面可以住一位学生。四个学生将会共用一个客厅，卫生间，淋浴等公共场所。笔者当年住在 North Tower 的 22 层，通过窗外就可以看到 UTown 的景色。从 PGP 搬到 UTown 的时候应该是 2014-07-15 这个日子，随后笔者就在 UTown 一直住到了博士毕业。不过当时从 PGP 搬过来的时候，论文中的核心部分已经解决，剩下的就只是修修补补的工作，只需要把论文写完即可。于是在博士生第五年的时候，除了日常去 S17 的 Office 之外，其余大部分时间都是在 UTown 度过的。相较于在 PGP 的时光，由于没有了课题的巨大压力，在 UTown 的时光整体来看还是十分舒适和惬意的。（未完待续）

Photo 14-7-18 10 55 26 Photo 14-7-18 10 56 06

National University of Singapore

NUS E-Open House

March 3, 2020 zr9558 Leave a comment

新加坡（Singapore）

新加坡被誉为“花园城市”，又称为“狮城”，它位于马来半岛的南端，地处在马六甲海峡最南端的位置。属于热带雨林气候的新加坡，一年只有雨季与旱季两个季节，从而导致狮城降雨频繁，气温长期处于 25-34 摄氏度左右。新加坡是一个多民族融合的国家，英语作为其官方语言。

新加坡国立大学（National University of Singapore）

位于新加坡西海岸的新加坡国立大学（National University of Singapore）是一所综合性大学，近十年在泰晤士报（Times）世界大学排名和 QS 世界大学排名上，NUS 都长期名列前茅，属于亚洲的顶级学府。NUS 的前身成立于 1905 年，今年恰好是 NUS 115 校庆的年份。发展至今，NUS 已经是拥有多个学院的综合性大学，其中包括理学院（Faculty of Science），工学院（Faculty of Engineering），商学院（Business School）等诸多学院。

NUS Open Day（新加坡国立大学开放日）

学校的排名除了全体教职工的努力之外，优秀的生源也是保障学校能够持续运营的必要条件。因此，在每个学年的第二个学期，大约在 2，3月份的时候，NUS 都会开展 Open Day，目的之一就是对新加坡的本地学生和国际学生宣传 NUS。除了宣传自己的硬件设施，专业特色，教授质量之外，还会宣传校园生活（Campus Life），基础设施和历年优秀学生的就业情况。

当年笔者还在 NUS 就读的时候，就有幸参加过 NUS 的 Open Day，不过当年只是为了凑热闹而去观看了一些活动，并没有在申请学校之前就享受到这些福利。每年到了 Open Day 之际，都会有大批中学生或者理工院校的学生前来 NUS，通过参观 NUS 的情况，来判断 NUS 是否适合自身的发展。

NUS E-Open House

由于 2020 年的开局实在是不利，世界上的诸多国家都受到了新型冠状病毒的影响，自然新加坡也不例外。可能是因为这个原因，2020 年的 NUS Open Day 就从线下（Offline）搬到了线上（Online），通过网络这一个重要的媒介来开展 Open Day。今年的 NUS E-Open House 应该是历史上首次在网上举办的校园开放日活动，同时校方在 NUS 的官网上也进行了大力的宣传。

NUS-E-OpenDay-Facebook-3 — NUS E-Open House 2020

除此之外，NUS 也通过 Facebook 账号也进行了推广，并且学生们可以通过 Facebook，YouTube，Instagram，Zoom 等诸多社交网络工具来全方位的了解 NUS。

由于本次的 NUS E-Open House 为期九天，从 2020-02-26 到 2020-03-05，个人感觉 Open Day 不如称之为 Open Week。在这段时间，NUS 的所有学院将会通过网络向全世界展示 NUS 的魅力。在这段时间里，学生和家长，包括想要继续深造的职场人士都可以通过社交网络了解到 NUS 的方方面面。

NUS-E-OpenDay-Web-1 — 2020 年 NUS E-Open House 的日程安排

除了学院之外，每个学院都会基于自身的条件来设计相关的项目。因此，向外部人士展示各个学院的项目细节也是一个非常重要的环节。

当然，学习是大学生活的一部分，除了学习之外，校园生活也是许多家长和学生关心的话题之一。于是，丰富多彩的校园生活则是 NUS 的一大特色。

NUS 的院系

Faculty of Arts and Social Sciences（艺术与社会科学学院）

前几年，随着直播（Live）这门技术的蓬勃发展，直播已经融入了人们生活的方方面面。FASS 的学生则是通过直播的方式向大家介绍了学院的地理位置，食堂，校园生活，学习等诸多内容，让观众感受到了在校学生的热情与活力。

NUS-E-OpenDay-FASS-5

当然，只有学生的宣传是远远不够的，此时需要整个院系的努力才能够把宣传力度加到最大。在 NUS e-Open House 的主页上，可以找到 FASS 的活动时间和宣传安排，使用各种各样的宣传渠道，对这个专业感兴趣的学生就会早早订阅并且关注这一消息。

在 FASS 中，学生们可以选择的专业也是非常多的，不仅包括中文，日语，哲学，历史等常见课程，还包括心理学，社会科学等方向。

Business School（商学院）

在商学院，不仅邀请了院系的教授与未来的学生们直接对话，还邀请了学生代表进行发言，通过亲身经历来向大家展示这几年在商学院得到的知识与心得。

NUS-E-OpenDay-Business-3 NUS-E-OpenDay-Business-1 NUS-E-OpenDay-Business-2

School of Continuing and Lifelong Education（持续与终身教育学院）

学习不仅仅是学生时期的事情，在这个高速发展的时代，学习这件事情将会伴随我们每个人一生的时间。于 2016 年成立的 School of Continuing and Lifelong Education 会给每个成年人提供持续教育的机会。它主要是为了想学习新技能的成年人或者想拿到学位的成年人而设置的。在这个学院有本科项目和硕士项目，还有短期项目或者各种培训。如果想拿到某些证书或者学位，其实这里是个不错的选择。

当年在新加坡国立大学（NUS）读书的时候，身边就有一些攻读 Master 学位的同学，当时他们就已经在公司工作了，在工作之余会选择 part time 的硕士在进行攻读。只要达到了学分或者做完了相应的 Project 就可以拿到 Master 学位。

NUS-E-OpenDay-SCALE-7

NUS-E-OpenDay-SCALE-2

Department of Mathematics（数学系）

新加坡国立大学数学系的前身可以追溯到 1929 年的 Raffles College。当时理学院开设了数学，化学，物理三门课程，不过总共也就只有十个学生和三位教师，其中有一位是数学教师。第一届数学系的领导（从1931年到1959年）是 Alexander Oppenheim 教授，他是在美国芝加哥大学获得的博士学位。从 1929 年开始，在新加坡的教育系统中，数学教育事业得到了巨大的发展，对现在的新加坡国立大学和南洋理工大学的建立起到了至关重要的作用。

随着 NUS 的建立，数学系就进入了一个新的时代。新的校区在 Kent Ridge，1986 年理学院和数学系就在这里成立。这个时候，数学系就有了巨大的发展，不仅在本科生的招生规模方面有了巨大提高，在研究生项目规模上也有了一定的深度的提升。

在本次的 NUS E-Open House 中，Science 的每个院系都提供了相应的宣传资料，并且所有学科统一模板，由各个学科的教授来向大家介绍这些专业的背景，学习内容和相关优势。在 NUS 的数学系（Department of Mathematics），则是由 Prof Tan 来给大家介绍院系的相关内容。

NUS-E-OpenDay-Science-3

通过 Prof Tan 的PPT，我们可以得到选择数学系的五大原因分别是：

Finding Good Jobs；
Wide Range of Career Choices；
Multiple Pathways；
Myriads of Real World Applications；
Life Long Learning Skills。

NUS-E-OpenDay-Science-Math-6 — 数学系学生的就业统计

从学校统计到的数学系学生毕业出路来看，其实数学系学生们的选择范围还是相对广泛的。不仅包括教育行业（Education）和科研行业（Research），还包括金融（Finance），科技（Technology），信息管理（Information Management）等诸多热门方向。就业的公司不仅包括星展银行（DBS），花旗银行（CitiBank），还包括 Google，Facebook 等科技公司。如果学生在读大学之前并不确定未来要做什么方向的话，并且对学习数学有一定的能力和兴趣的前提下，在本科期间选择数学专业其实是一个不错的选择。

NUS-E-OpenDay-Science-Math-8 — 数学系学生的就业方向

俗话说，数学是科学的基础。数学不仅仅是数学书上的一道道定理，而是可以解决现实生活问题的重要工具。在天气预测方向，动力系统（Dynamical System）有着独特的应用；在机器学习领域，微积分，线性代数与概率论则为这门学科提供了理论基础。在金融领域，Black-Scholes 方程，Monte Carlo 模拟则是其中的重要模型。从这些学科的发展来看，数学不仅可以为这些学科提供理论基础，也是解决这些学科难题的重要工具之一。

NUS-E-OpenDay-Science-Math-12 — 数学在各个学科中的应用

NUS-E-OpenDay-Science-Math-16 — NUS 数学系的项目

数学系的课程以难著称，有人说：“没有比数学书更好的劝退材料”。要想劝退一个人学习数学，其实是非常简单的。在 NUS 数学系，数学系的课程难度就是课程 ID 的第一个数字。例如 MA1100，就表示 Level 1 的课程，属于数学系的入门课。MA4270 就表示 Level 4 的课程，属于数学系本科课程中较难的课程。学生们可以根据自身的情况和时间安排选择相关的课程进行学习。

当然，每个专业都有着不同的细分方向，虽然刚入学的时候，同一院系的学生所学课程都是一样的，但是随着时间的迁移，不同的学生会选择不同的方向进一步的深造。在数学系，通常会有 Math（数学），Applied Mathematics（应用数学），Quantitative Finance（金融数学），数据科学与分析（Data Science & Analytics）这几个方向。学生们需要在学习的过程中根据自身的兴趣和需要选择最合适自己的方向。

NUS-E-OpenDay-Science-Math-17 — NUS 数学系的项目

在本科的时候，通常都会有不少学有余力的学生，不满足于只学习一门课程，或者不满足于仅仅拿到一个学位。此时，NUS 数学系还可以提供双学位的项目（Double Degree Programme），或者选择双主修的项目（Double Major Programmes）。学生可以在其中选择经济学，管理学，商业分析，计算机科学，信息安全这几个方向进行辅修，从而把在数学系学到的知识进一步地应用在其他学科中。

Application（NUS 的申请）

在新加坡国立大学的官网上，可以找到其申请的入口，学生们可以根据自身的情况申请相应的项目，从而度过一个丰富多彩的大学生活。

结束语

通过新加坡国立大学的校园开放日（NUS Open Day），我们可以了解到 NUS 的院系特点，校园生活，学习项目等很多内容。而 NUS E-Open House 则是在这个特殊时期的一个创新，充分利用互联网的优势，把 NUS 的特色通过网络传递给大家，向大家充分展示了 NUS 的魅力之处。

参考资料：

职场纵横, 每天的时间。。。。

为什么【在家办公比上班还累】？

February 8, 2020 zr9558 Leave a comment

野原广志在家工作的一天

《在家里工作的爸爸哦》

看到这里想必大家也知道为什么在家里工作效率低了。从《蜡笔小新》这个纪实版动画可以分析得知大致有以下几个原因：

工作与生活的时间界限感不强，导致工作的开始时间不明确；
被家里的琐事频繁打断：例如洗衣服，照顾小孩，收快递，送东西等；
自制力不强：容易被电视，电脑，游戏等好玩的东西所吸引。

因此，野原广志在家一天工作的产出几乎是零。

解决这类问题的思路有以下几个：

1. 划清工作和生活的界限：找一个单独的房间，不被打扰，专心致志地工作。让自己从工作切换到生活的成本增加，而不是处于可以随意地切换工作与生活的状态。如果找不到相应的房间，就找一个桌子，日常生活不要坐到桌子前。但是一旦坐到了桌子前，就表示进入了工作状态。就表示不能够轻易地被打扰，要尽快开始一天的工作。

2. 使用日历表：根据自身情况，制定基于自己工作和生活状态的日历表，使用 Google Calendar 或者其他计时工具，把自己一天的日程安排下来，让自己在固定的时间内做固定的事情，不要被其他事情所打扰。

例如下图，就是笔者在上班的时候的日程安排：

3. 制定详细的工作计划：包括制定每天需要输出的工作内容，无论是写了多少企划书，写了多少页 PPT，写了多少代码，完成了多少个功能点，都需要有非常明确的指标。只有制定了明确的目标，才能够督促自己完成任务。

4. 规划固定时间：吃饭，睡觉，娱乐等固定的时间先划分好，不到万不得已，不能随便改变生活规律，一旦改变生活规律，要想纠正回正常状态就会非常的困难。

5. 区分工作内容：提前区分哪些是需要团队协作的项目，哪些是靠自己就能够独立完成的项目。对于需要团队协作的项目，尽量提前安排好时间，毕竟跨地域沟通需要一定的成本，需要对方也在线，把时间和人提前安排好会明显优于临时安排。对于那些靠自己就能够独立完成的项目而言，只需要保证工作的时候不被打扰，专心致志地做完即可。

职场人士在家工作确实可以体验自由职业者谋生或者博士生做科研的状态，确实也可以做到自行安排自己的时间，但是要做时间的主人，而不是做时间的奴隶。

Number Theory

素数之美

February 2, 2020 zr9558 Leave a comment

素数的定义

素数是在中小学课本里面就会出现的数学概念，它指的是只能够被 1 和它本身整除的正整数。在正整数中，2, 3, 5, 7, 11 等都是素数。同时，每一个正整数（不小于 2）都可以写成多个素数的乘积，例如 $35 = 5 * 7, 10 = 2 * 5$ 。从素数的定义可以看出，判断一个数是否是素数是需要通过“乘法”的。而在数学的研究历程中，数学家们同样也关心由素数之间的加法所产生的奇妙结论。

哥德巴赫猜想（Goldbach’s Conjecture）

随着徐迟的报告文学《哥德巴赫猜想》的问世，哥德巴赫猜想在国内早已家喻户晓。其中，哥德巴赫猜想包括两个部分：

[Theorem] 每一个大于 7 的奇数都可以写成三个素数之和；
[Conjecture] 每一个大于 6 的偶数都可以写成两个素数之和。

从猜想的陈述来看，如果第 2 部分是正确的，那么可以根据公式 $n = (n-3) +3$ 直接得到第 1 部分是正确的，因此第 2 部分被称为强哥德巴赫猜想，第 1 部分被称为弱哥德巴赫猜想。其中哥德巴赫猜想的第 1 部分已经被彻底解决，而哥德巴赫猜想的第 2 部分目前最好的结果被称为陈氏定理（ Chen’s Theorem）。用数学的语言来说，这两个定理的陈述分别是：

[Theorem (Vinogradov)] 假设 $N$ 是一个奇数，令 $r(N) = \sum_{p_{1}+p_{2}+p_{3}=N}1$ 表示关于 $N$ 的计数函数，其中 $p_{1}, p_{2}, p_{3}$ 都是素数。则存在一个一致有界的函数 $\Omega(N) \in (c_{1},c_{2})$ （ $c_{2}>c_{1}>0$ ）对于充分大的奇数 $N$ ，有以下式子成立

$r(N) = \Omega(N)\cdot \frac{N^{2}}{(\ln(N))^{3}}\cdot\bigg\{1+O\bigg(\frac{\ln\ln(N)}{\ln(N)}\bigg)\bigg\}.$

备注：从以上公式可以看出， $\lim_{N\rightarrow \infty} r(N) = +\infty.$ 换句话说， $\exists N_{0}, \forall N\geq N_{0},$ 弱哥德巴赫猜想成立。

[Theorem (Chen)] 假设 $N$ 是一个偶数，令 $r(N)=\sum_{p+n=N}1$ 表示关于 $N$ 的计数函数，其中 $p$ 是素数， $n$ 表示最多为两个素数的乘积。则当 $n$ 充分大的时候，有以下式子成立：

$r(N) >> \Omega(N)\cdot \frac{2N}{(\ln(N))^{2}},$

其中 $\Omega(N) = \prod_{p>2} \bigg(1-\frac{1}{(p-1)^{2}}\bigg)\prod_{p|N, p>2}\frac{p-1}{p-2}.$

备注：

在哥德巴赫猜想的研究过程中，通常数学家把偶数可表示为 $a$ 个素数的乘积与 $b$ 个素数的乘积之和这个问题，简称为 $a + b$ 问题。所以，陈景润证明的 “1+2” 并不是指 1+2 = 3，而指的是对于每一个充分大的偶数，要么是两个素数之和，要么是一个素数加上两个素数之积。其实可以简单的理解为 $p_{1}+p_{2}$ 或者 $p_{1}+p_{2}\cdot p_{3}$ ，在这里 $p_{1},p_{2},p_{3}$ 都是素数。从以上公式可以看出， $\lim_{N\rightarrow \infty} r(N) = +\infty.$
1920 年，挪威数学家 V.Brun 证明了 “9+9″，开启了数学家研究哥德巴赫猜想之路；1966 年，中国数学家陈景润证明了 “1+2″，把素数的筛法推向了顶峰。

孪生素数猜想（Twin Primes Conjecture）

在上千年的素数研究历程中，除了哥德巴赫猜想，孪生素数（Twin Primes）的研究也是数论中的一个重要课题。所谓孪生素数就是相差为 2 的两个素数，例如 $(3,5), (5,7), (11,13),\cdots$ 等等。因此，就有人提出猜想：孪生素数有无穷多对。换句话说，如果用 $p_{n}$ 表示第 $n$ 个素数，那么孪生素数猜想就是 $\liminf_{n\rightarrow +\infty}(p_{n+1}-p_{n})=2$ . 除了孪生素数本身之外，也有学者猜测，对于所有的正整数 $k\geq 1,$ 形如 $(p,p+2k)$ 的素数对同样有无穷多对。于是，在网上就有人对于有限的素数对进行了计算，让大家更好地看到素数之间的分布情况。

下面是部分关于素数间距（小间距，Small Gaps）的结论：

1940 年，Paul Erdos 证明 $\exists c>0$ 使得 $\liminf_{n\rightarrow\infty} \frac{p_{n+1}-p_{n}}{\ln(p_{n})}<c.$
2005 年，Daniel Goldston，Janos Pintz 和 Cem Yildirim 证明 $\liminf_{n\rightarrow\infty}\frac{p_{n+1}-p_{n}}{\ln(p_{n})}=0.$
2007 年，上述结果被改进为 $\liminf_{n\rightarrow\infty}\frac{p_{n+1}-p_{n}}{\sqrt{\ln(p_{n})}\cdot (\ln\ln(p_{n}))^{2}}=0.$
2013 年，张益唐证明了 $\liminf_{n\rightarrow\infty}(p_{n+1}-p_{n}) < 7 * 10^{7}$ ，随后这个结果被改进到 246。

除了素数之间的小间距之外，素数之间的大间距（Big Gaps）同样也有很多结论：

1931 年，Erik Westzynthius 证明 $\limsup_{n\rightarrow\infty}\frac{p_{n+1}-p_{n}}{\ln(p_{n})} =\infty.$
2014 年，Kevin Ford, Ben Green, Sergei Konyagin, Terence Tao 和 James Maynard 证明 $p_{n+1}-p_{n}>c\cdot \frac{\ln(n)\cdot \ln\ln(n) \cdot \ln\ln\ln\ln(n)}{\ln\ln\ln(n)}$ 对于某个 $c>0$ 和无穷个 $n$ 成立。

素数定理

在研究素数的过程中，研究素数的分布规律就是这一切的关键所在。其中，素数定理则是描述素数分布的一个重要结论。类似的，关于孪生素数的分布也有一个上界的估计。

[素数定理] 假设 $\pi(x)$ 表示不大于 $x$ 的所有素数的个数，那么 $\lim_{x\rightarrow \infty}\pi(x)/(x/\ln(x)) = 1.$

[孪生素数个数的上界] 假设 $\pi_{2}(x)$ 表示不大于 $x$ 的所有孪生素数个数，那么存在常数 $C>0$ 使得 $\pi_{2}(x)\leq C\cdot x/(\ln(x))^{2}.$

备注：从这两个定理可以粗糙地刻画出素数与孪生素数在实数轴的分布情况，并且可以看出孪生素数相对于素数则是少很多的。因为 $\lim_{x\rightarrow+\infty}\pi_{2}(x)/\pi(x)=0.$

素数的性质

在中小学的竞赛部分，大家总能够接触到一个关于素数的定理。

[Theorem (Euclid)] 素数有无穷多个。

证明：假设素数是有限个，不妨设为 $p_{1},\cdots, p_{n}$ ，那么 $N = \prod_{i=1}^{n}p_{i}+1$ 就是合数，但是它却不能被所有的素数 $p_{1},\cdots,p_{n}$ 整除，所以导致矛盾。因此素数是无穷多个。证明完毕。

除此之外，在大学里面学习级数的时候，通常都会研究调和级数（Harmonic Series）的性质。所谓调和级数指的就是所有正整数的倒数和，形如：

$S(x) = \sum_{1\leq n\leq x} \frac{1}{n}.$

从定积分与级数的关系可以得到 $\lim_{x\rightarrow +\infty}S(x) = +\infty$ 并且 $S(x) = \ln(x) + O(1).$ 也就是说，所有正整数的倒数和是发散的。

利用这种思路，其实可以分析所有素数的倒数和，也就是说 $\sum_{p \text{ prime}} \frac{1}{p}.$ 通过欧拉公式可以得到：

$\sum_{n=1}^{\infty}\frac{1}{n} = \prod_{p\text{ prime}}\bigg(1+\frac{1}{p}+\frac{1}{p^{2}}+\cdots\bigg)= \prod_{p\text{ prime}} \frac{1}{1-\frac{1}{p}},$

两边取对数可以得到 $\ln\bigg(\sum_{n=1}^{\infty}\frac{1}{n}\bigg) = \sum_{p\text{ prime}} -\ln\bigg(1-\frac{1}{p}\bigg),$

由于 $-\ln(1-x) = x + O(x^{2})$ ，并且 $\sum_{n=1}^{\infty}\frac{1}{n}=+\infty$ , $\sum_{n=1}^{\infty}\frac{1}{n^{2}}=\frac{\pi^{2}}{6},$ 可以得到

$\ln\bigg(\sum_{n=1}^{\infty}\frac{1}{n} \bigg)= \sum_{p\text{ prime}}\frac{1}{p} + O\bigg(\sum_{p\text{ prime}}\frac{1}{p^{2}}\bigg).$

等式的左边是发散的，右侧的第二项是收敛的，因此右侧的第一项（素数的倒数和）是发散的。进一步地，可以得到两个结论：

$\sum_{p\text{ prime}, p\leq x} \frac{1}{p} = \ln\ln(x) + O(1);$
$\prod_{p\text{ prime}, p\leq x}\bigg(1-\frac{1}{p}\bigg)^{-1} = c\cdot\ln(x)+o(1),$ 这里， $c>0$ 是一个常数。

至此，我们得到了两个级数的定理：

[Theorem] 所有正整数的倒数和是发散的；
[Theorem] 所有素数的倒数和是发散的。

从第 2 个结论同样可以得到素数是无穷多个。于是，就有数学家猜测如果孪生素数的倒数和是发散的，那么孪生素数同样也是无穷多对。但是在 1915 年，数学家 Brun 证明了，孪生素数的倒数和是收敛的，这个收敛的数字也被称为 Brun 常数。

[Theorem] 所有孪生素数的倒数和是收敛的。

证明：通过孪生素数个数的上界公式，可以得到存在 $C>0$ 使得对于充分大的 $x$ ，有

$\pi_{2}(x) \leq C\cdot \frac{x}{(\ln(x))^{2}}$

成立。假设素数序列 $p'_{1}, p_{2}',\cdots, p_{n}',\cdots$ 使得 $p, p+2$ 都是素数，那么 $n = \pi_{2}(p_{n}') \leq C\frac{p_{n}'}{(\ln(p_{n}'))^{2}}\leq C\frac{p_{n}'}{(\ln(n))^{2}}$ ，进一步可以得到

$\frac{1}{p_{n}'} \leq C\frac{1}{n\cdot (\ln(n))^{2}}$

对于充分大的 $n$ 成立。而右侧是收敛的，i.e. $\sum_{n=2}^{\infty} \frac{1}{n\cdot(\ln(n))^{2}}<+\infty.$ 因此，孪生素数的倒数和是收敛的。证明完毕。

备注：由于孪生素数的倒数和是收敛的，因此，通过孪生素数的倒数和来证明孪生素数有无穷多对这条路就被封死了。

在研究孪生素数的过程中，其目的是为了研究素数之间的间距究竟能有多小，也就是分析 $\liminf_{n\rightarrow+\infty}(p_{n+1}-p_{n})$ 的上界。同样的，也可以研究素数之间的间距究竟有多大，并且可以分析其量级大约是多少，此时就需要研究 $\limsup_{n\rightarrow+\infty}(p_{n+1}-p_{n}).$

[Theorem] 对于充分大的 $x$ 而言，在 $[1,x]$ 内，素数之间的最小间隔 $\min_{p_{n}\leq x} (p_{n}-p_{n-1})\leq (1+o(1))\ln(x);$ 同时，素数之间的最大间隔 $\max_{p_{n}\leq x}(p_{n}-p_{n-1})\geq (1+o(1))\ln(x).$

证明：考虑区间 $[1,x]$ ，通过素数定理可以得到在 $[1,x]$ 区间内的素数大约是 $\pi(x) \sim x/\ln(x)$ 个。于是把该区间 $[1,x]$ 切割成长度为 $\ln(x)$ 的子区间，区间的个数为 $x/\ln(x),$ 通过鸽笼原理 (Pigeonhole Principle) 可以得到此定理的结论。

备注：除此之外，证明相邻素数的间隔没有上限还可以用构造法。考虑 $n!+2, n!+3, \cdots, n!+n$ 这 $n$ 个连续的合数，所以两个相邻的素数必在 $[n!+2, n!+n]$ 这个区间两侧。因此相邻素数的间隔没有上限，i.e. $\limsup_{n\rightarrow+\infty}(p_{n+1}-p_{n})=+\infty.$

Eratosthenes 筛法（Eratosthenes Sieve Method)

Eratosthenes 筛法是数学家 Eratosthenes 提出的一种筛选素数的方法，其思路比较简单：想要筛选出 $[2,n]$ 中的所有素数，则首先把 $[2,n]$ 中的所有正整数按照从小到大的顺序 $2, \cdots, n$ 来排列，然后按照如下步骤执行：

读取数列中当前最小的数 2，然后把 2 的倍数全部删除；
读取数列中当前最小的数 3，然后把 3 的倍数全部删除；
读取数列中当前最小的数 5，然后把 5 的倍数全部删除；（4 已经被第一步去掉了）
读取数列中当前最小的数 7，然后把 7 的倍数全部删除；（6 已经被第一步去掉了）
循环以上步骤直到 $[2,n]$ 中所有的数被读取或者被删除。

其算法复杂度为 $O(n\ln(n))$ 。

Brun 筛法（Brun Sieve Method）

在数学界发展出各种筛法，其重要目的之一就是为了解决孪生素数猜想和哥德巴赫猜想。除了 Eratosthenes 筛法之外，数学家 V. Brun 也发现了一种筛法，后人称之为 Brun’s Sieve。其目的就是为了估计孪生素数的上界，进一步得到计算孪生素数的倒数和。其主要结论就是 $\pi_{2}(x)\leq C\frac{x}{(\ln(x))^{2}},$ 其中 $C>0$ 是一个常数，并且 Brun 通过其筛法可以得到哥德巴赫猜想中的 “9+9″，在哥德巴赫猜想的发展中属于里程碑式的工作。

Question. 研究素数究竟有什么用？

Answer. 为了人类智慧的荣耀。

参考文献：

Small and Large Gaps Between Primes, Terence Tao, Latinos in the Mathematical Sciences Conference, 2015.
Bounded Gaps Beween Primes, Yitang Zhang, 2013.
Additive Number Theory, Melvyn B.Nathanson, GTM 164.
http://mathworld.wolfram.com/TwinPrimes.html.

PHD的生涯

在新加坡的这五年—生活篇（四）

January 29, 2020 zr9558 Leave a comment

在之前的文章里，笔者分别对新加坡留学生的日常生活，人物轶事，过年感受做了回顾。本次将会提到一个重要的话题，那就是电子游戏（Electronic Games）。

对于基础学科的博士生来说，在第一年都会修各种各样的课程，通常来说毕业之前需要修 8 门课，再加上一门听讲座的课程。并且在前两年参加博士生资格考试（Qualify Exam，分为 Writing 和 Oral 两个部分）。一旦博士生通过了博士资格考试，就开始全面进入科研阶段；否则就只能够拿硕士学位离开。众所周知，读书是很多人一起参与的事情，通过考试是可以分出高下的。唯独科研这一件事情，几乎就是自己一个人的战斗。在这一个人的漫漫科研长路中，虽然大部分时间笔者也在奋力向前，努力毕业，但偶尔难免心生疲惫。在疲惫或者想放松的时候，就自然而然地想到了电子游戏（Electronic Games）。

在中学的时候，笔者接触过《帝国时代》（Age of Empires），《暗黑破坏神2》（Diabol 2）等系列游戏；在大学期间，则开始玩《魔兽争霸3》（WarCraft 3），《Dota》；到了博士生的时候，接触到的游戏就更多了，包括但不限于《暗黑破坏神3》（Diabol 3），《炉石传说》（HeartStone），《坦克世界》（World of Tanks），《Fifa 2015》。

在中学的时候，虽然也接触过《暗黑破坏神 2》，但是由于高考等种种因素，没有认认真真地玩过这款游戏。真正开始玩这款游戏应该是在大四下学期（ 2009 年），当时 05级数学系的所有学生都已经从浦口校区搬到了鼓楼校区，男生们都住在筒子楼里面。之所以是筒子楼，意思就是学校的宿舍不够用了，拿家属区的房子来充当宿舍。其实在 2009 年之前，当时还在浦口的时候，暴雪公司就早已对外公布了《暗黑破坏神 3》的一些画面和角色，包括野蛮人的装备配置，闯关场景等诸多画面。但直到 2012 年 5 月份，暴雪公司才正式对外发售这款游戏。在对外发售这款游戏之前，NUS 数学系的很多博士生都沸腾了，纷纷通过各种渠道购买这款游戏，有的人提前在实体店预定，有的人提前在网上购买。笔者当时在暴雪公司的官网上花了 90 SGD 购买这款游戏，包括亚服，美服，欧服三个服务区。虽然博士生的学习十分繁忙，但是 NUS 的暑假总是来得特别早，每个学年基本上 4 月份就完成了教学工作，5 月份的话就已经开始了暑假的生活，直到 7 月底才会结束假期，开始新学年的生活。在种种因素下，除了必要的工作和科研之外，有的博士生就开启了打暗黑 3 的生活。

玩 PC 游戏比较依赖电脑的性能，当时虽然花钱购买了游戏，不过当时的电脑是 2008 年的电脑，在鼓楼珠江路购买的 Lenovo 的 IdeaPad Y 530。即使此电脑性能不行，但也满足了《暗黑破坏神 3》这款游戏的最低配置。因此，在身边的同学都纷纷入坑的情况下，笔者也跟着大家开始玩这款游戏。

LenovoidealpadY530 — Lenovo IdeaPad Y530

不过到了 2010 年之后，iPhone 开始风靡全球，智能手机开始逐渐走入大家的视野。到了 2011 年 7 月份的时候，笔者决定买一台 iPhone。在新加坡，如果是本地人或者 PR，或者有工作签证的人，购买一台 iPhone 的合约机是非常简单的。每人可以签约一台手机，预购费貌似是 200 SGD 左右，然后每个月大概也就花 50 SGD 左右就可以将其拿到。但是对于持有学生签证的人而言就没有那么方便了，必须要缴纳 800 SGD 的押金，时长两年，才能够把手机拿到手。不过在当时 iPhone 的热潮下，跟着贾博士，尹博士一起去 vivo city 的 StarHub 办理手续。除了钱之外，还需要当时银行来往的信件，因为信件上面有每个人的通讯地址，每个月 StarHub 会通过这个地址给签约人寄送当月的通讯费用。由于没有带银行来往的信件，三人还往返了一趟学校与 vivo city，跑了两次才将手机拿到手。

“一入苹果深似海，从此钱财是路人”。苹果的产品确实做得很好，除了 iPhone 4 之外，后续笔者每年都会入手一些苹果的其他产品，包括 iPod Shuffle，itouch，键盘，鼠标，iPad，Mac Air。除了一体机 iMac 没有入手之外，其他的苹果产品均已体验过。当时基于这些产品，在 2013 年的时候把自己的电脑从 Lenovo 换成了 Mac Air，并且在博士生的办公室（Graduate Office 5）搭建了自己的工作环境（自然也可以无缝切换成游戏环境）。

当年在新加坡的时候，师门总共有四个人，分别是大师兄，二师兄，三师兄再加上笔者。其中，大师兄从来不玩任何游戏，二师兄与三师兄在那段时间都会玩《暗黑破坏神 3》。特别地，三师兄是众多博士生中唯一一个能够通过游戏本身回本的。在刚推出的时候，的《暗黑破坏神 3》有拍卖场，很多玩家都会在上面拍卖自己获得的装备，武器，金钱等道具。三师兄又特别擅长研究各种各样的游戏攻略和刷顶级道具，利用拍卖场这一套现工具，三师兄在玩了几个月之后顺利地把自己的装备套现，并且卖到了和游戏本身差不多的价钱。除了博士生同学之外，当年笔者还住在 West Coast Plaza 的 Block 602 #11-28，房东就是 NUS 的讲师，也一起入手了暗黑 3 这款游戏。可见暗黑系列游戏在 80 后一代人心中的地位。

随着时间的推移，科研任务的逐渐加重，新学期的重新开始，博士生中的暗黑 3 热潮已经逐渐淡去。到了 2014 年，暴雪的另一款游戏引起了我的注意，那就是《炉石传说》。炉石传说与暴雪的其他策略竞技类游戏不一样，它是一款卡牌游戏。根据卡牌的抽取顺序和英雄技能的不同，有着完全不一样的打法和策略。刚入门的时候门槛不高，笔者也能够有一些胜率，但是到了某天晚上，不知道是状态不好还是水平不够，一直在输，总没有克制对手的方法。于是，当天晚上笔者就往炉石传说里面充了钱，购买了很多卡包，当晚就拥有了很多卡牌，后来胜率又逐渐提升。果然，游戏只要充钱就会变得更强。

到了 2014 年的下半年，笔者也忙于撰写毕业论文和教学任务。不少当年在一起玩游戏的学长们也已经毕业，即使没有毕业也在努力搞论文。此时也很难凑到一批人像当年一样打暗黑 3。同时，随着智能设备的兴起，电子游戏也已经从以前的 PC 端转到了移动端，包括 iPad 和手机端。在手机端上，其实是没有办法进行复杂的游戏操作的，暗黑 3 这种游戏在 iPad 和手机上几乎是不可行的。在撰写毕业论文的同时，为了舒缓压力，就又挖掘了两款新的游戏，分别是《坦克世界》（World of Tanks）和《Fifa 2015》。

在坦克世界里，相较于各种坦克，笔者更倾向于用反坦克炮，靠灌木丛或者遮掩物躲在一处。一旦队友发现了敌方的踪迹就瞄准，然后就用反坦克炮给予一记重炮。其中，虽然苏系坦克的移动速度一般，但是火力生猛，深受笔者喜爱。

在本科的时候，笔者就经常和同学们一起玩实况足球，当时流传的一句话就是“无兄弟，不实况”，意思就是实况足球这种游戏确实要与同学一起竞技才能够体现其游戏乐趣。否则玩足球游戏就是以通关为目的，只要赢了或者获得冠军，游戏就已经结束了。当时大约是 2014 年的 12 月，2010 级的博士生要么已经提交了论文，要么正在提交论文的路上，笔者很少有机会与大家玩网络游戏，当时玩得更多的是单机版游戏，或者就与网上不认识的人随便打一下。

在读博士这几年，关于电子竞技或者生活的回忆还是很多的。生活嘛，总要找一点乐趣才能够坚持下去（未完待续）。

PHD的生涯

在新加坡的这五年—生活篇（三）

January 28, 2020 zr9558

在上一篇文章《在新加坡的这五年—生活篇（二）》中，笔者回顾了那些年在新加坡读书的时候，所遇到的一些有趣的人和事。本篇文章会回顾之前在新加坡的过年轶事。

最近几天正值新春佳节，虽然 2020 年的开局实在是不够顺利，无论是肆意传播的疫情，还是科比的突然离世，都给新的一年带来了不少沉重的色彩。但是新年总要有新的气象，每年都要立一些新的 Flag（虽然 2010 年的 Flag 可能至今都没有完成）。

li_flag

如果是在国内过春节，除了一些特殊的原因之外，那自然是需要回家乡与亲人团聚，坐在一起吃年夜饭，共同看春节联欢晚会。但是对于留学生而言，回国与家人共度春节就显得没有那么方便了。虽然新加坡的工作人士有很多的年假，但是他们的公共假日却少得可怜。很多外国人也很难想象中国人可以有春节黄金周，十一黄金周这种一整周的假期。在新加坡的官网上可以看到，新加坡的公共假日基本上都是一两天，除了中国农历新年放假两天之外，元旦新年，劳动节，新加坡国庆日，圣诞节等诸多假日都只有一天的假期。在这些公共假日里，从宿舍区去办公室是没有校车的，只能靠自己的 11 路走来走去。

SingaporePublicHoliday2020 — 2020 新加坡的公共假日

新加坡是一个多民族融合而成的国家，虽然农历新年对于中国人来说意义非凡，但是私底下有的老师对于留学生春节期间离开新加坡返回中国过节这件事情是有一定的意见的。毕竟这段时间学校已处于开学阶段，一旦有学生要离开助教（Teaching Assistant）或者助研（Research Assistant）岗位，该学生或者院系就需要找到另外的学生来临时代替他的工作，确实会对教学和科研造成一定的影响。不过这个只是说不建议留学生回国，其实留在新加坡的话留学生们也不需要来办公室，安心在家休假过节即可。对于很多留学生而言，即使不能够回国与家人欢度除夕新年，在新加坡体验异国他乡的春节其实也还算一个不错的选择。新加坡是以华人为主的社会，农历新年在新加坡绝对算得上是一个非常隆重的节日。每到农历新年之际，无论是 Chinatown，圣淘沙，还是每家每户，都是张灯结彩，充满了喜庆的氛围。

从 NUS 每年的学校日历来看，对于身处新加坡的留学生而言，农历新年的时候正好是刚刚开启新的学期的时候。即使回国也只能迅速返回学校，毕竟只有两天的假期。除了少数回国的学生之外，大部分留学生们都会在新加坡自行组织过农历新年的活动。对于学校而言，虽然是正月初一和正月初二放假，其实到了大年三十下午的时候，就已经不再组织教学活动了，因为就算上课也会有不少学生缺席。无论是学校的教职工，还是学生，都会在除夕当天下午准备过年的食物和喜庆的礼品。在正月初一和正月初二这两个特殊的日子里，很多商店，超市，饭店都是不开门的，即使是每天营业时间最长的昇菘超市（shengsiong），也会选择在大年三十的下午 3:00 提前下班，让员工们提前回去过除夕。因此，提前去超市准备过春节的食物和年夜饭就是所有留学生的必做的事情。

NUS_2019_2020_calendar — NUS 的 2019-2020 年度日历

要想从 NUS 的 UTown 或者 PGPR 两个宿舍区去昇菘超市可没有那么容易，除了等待 West Coast Plaza 的穿梭巴士和 NUS Bus，只能够选择步行的方式，从数学系 Block S17 出发，步行穿过工学院（Faculty of Engineering），再翻山穿过 Clementi Woods，最后走到昇菘超市。去时容易，返回困难。去程的时候双手空空荡荡，返程的时候可是大包小包全部装满了食物，饮料和酒水。通常来说，这往返一趟办公室和昇菘超市，需要两个小时左右的时间。不过，作为年夜饭的准备，还是非常值得的。

从 2011 年算起到 2014 年，笔者前前后后在新加坡度过了四个春节，那几年的除夕夜分别是在 West Coast Plaza 的川江号子，UTown，UTown，ChinaTown 的国府珍锅吃的年夜饭。印象最深刻的就是在 2014 年大年三十的时候，吃完年夜饭，回到 PGP 之后打开电脑，写下了这一段话，立下了 2014 年的 Flag：

今年是农历的除夕，也就是大年三十，依旧和PHD们出去吃了顿火锅。不知道是不是最后一次在新加坡过年了。前前后后也过去四年了。第一次是West Coast Plaza的川江号子，第二次是utown宿舍，第三次也是在utown宿舍，第四次是在Chinatown的国府珍锅。不知道明年还有没有机会在新加坡过年了。最希望的就是尼玛交了论文就直接可以考虑毕业了，就不用在新加坡这个地方过年了。

Jan 31, 2014 00:18

于是，为了兑现 2014 年初的 Flag，笔者在 2014 年整整奋发图强了半年，终于在 2014 年 7 月份的时候把论文的主要步骤自行推导出来，在 2014 年下半年开学之后花了两三个月时间把所有步骤整理并且装订成册，在 2015 年 1 月份的时候就在系统上提交了自己的博士论文。

每次提到毕业的事情都会觉得心情沮丧，所幸的是最终还是顺利地拿到博士学位。回想起博士生活的点点滴滴，除了科研的过程比较苦逼之外，其他的各种生活还是丰富多彩的。（未完待续）

Mathematics

拉格朗日四平方和定理

January 26, 2020 zr9558 Leave a comment

拉格朗日四平方和定理

每个正整数均可表示为四个整数的平方和。

Every positive integer is the sum of four squares.

例如：

$1=1^{2}+0^{2}+0^{2}+0^{2}$ ，
$2 = 1^{2}+1^{2}+0^{2}+0^{2}$ ，
$7 = 2^{2}+1^{2}+1^{2}+1^{2}$ 。

证明：可以直接验证如下恒等式

$(x_{1}^{2}+x_{2}^{2}+x_{3}^{2}+x_{4}^{2})\cdot(y_{1}^{2}+y_{2}^{2}+y_{3}^{2}+y_{4}^{2}) = z_{1}^{2}+z_{2}^{2}+z_{3}^{2}+z_{4}^{2}$ ，其中

$\begin{cases} z_{1}=x_{1}y_{1}+x_{2}y_{2}+x_{3}y_{3}+x_{4}y_{4} \\ z_{2}=x_{1}y_{2}-x_{2}y_{1}-x_{3}y_{4}+x_{4}y_{3} \\ z_{3}=x_{1}y_{3}-x_{3}y_{1}+x_{2}y_{4}-x_{4}y_{2} \\ z_{4}=x_{1}y_{4}-x_{4}y_{1}-x_{2}y_{3}+x_{3}y_{2}\end{cases}$

由于 1 与 2 都明显满足这个定理，那么只需要考虑大于 2 的正整数。而这些正整数都可以分解成素数的乘积，因此，只需要证明该定理对所有的素数成立，则使用以上恒等式就可以得到最终的结论。假设 $p$ 是一个奇素数。

由于 $\{a^{2}:a\in\{0,1,\cdots,(p-1)/2\}\}$ 里面有 $(p+1)/2$ 个不同的同余类， $\{-b^{2}-1: b\in \{0,1,\cdots,(p-1)/2\}\}$ 里面也有 $(p+1)/2$ 个不同的同余类，但是素数 $p$ 的同余类只有 $p$ 个，因此存在正整数 $a,b\in \{0,1,\cdots, (p-1)/2\}$ 满足 $a^{2}\equiv -b^{2}-1 (\mod p)$ 。也就是说 $a^{2}+b^{2}+1^{2}+0^{2}\equiv 0(\mod p)$ 。令 $n\in\mathbb{Z}$ 满足 $np=a^{2}+b^{2}+1$ ，则有 $p\leq np\leq 2(p-1)^{2}/4+1<p^{2}$ 。于是， $1\leq n<p$ 。

因此存在一个 $1\leq n<p$ 使得 $np = a^{2}+b^{2}+1^{2}+0^{2}$ 是四个整数的平方和。于是必定存在一个最小的正整数 $m$ 使得 $1\leq m\leq n<p$ 使得 $mp$ 为四个整数的平方和，不妨设为 $mp=x_{1}^{2}+x_{2}^{2}+x_{3}^{2}+x_{4}^{2}$ 。

Claim. $m=1$ 。

proof of the claim. 反证法，假设 $1<m\leq n<p$ 成立。令 $y_{i}=x_{i}(\mod m)$ 对于 $i\in\{1,2,3,4\}$ 成立，并且 $-m/2<y_{i}\leq m/2$ 。因此， $y_{1}^{2}+y_{2}^{2}+y_{3}^{2}+y_{4}^{2}\equiv(x_{1}^{2}+x_{2}^{2}+x_{3}^{2}+x_{4}^{2})\equiv mp \equiv 0(\mod m)$ 。令 $mr = y_{1}^{2}+y_{2}^{2}+y_{3}^{2}+y_{4}^{2}$ 。因此， $mr\leq 4(m/2)^{2}=m^{2}$ 。

如果 $r =m$ ，通过以上不等式得知 $r=m$ 等价于 $y_{i}=m/2$ 对于 $i\in\{1,2,3,4\}$ 都成立。此时， $mp = x_{1}^{2}+x_{2}^{2}+x_{3}^{2}+x_{4}^{2}\equiv 4(m/2)^{2} \equiv 0 (\mod m^{2})$ 。因此， $p$ 是 $m$ 的倍数，这与 $p$ 是素数， $m>1$ 矛盾。所以， $r<m$ 成立。i.e. $1\leq r<m\leq n<p$ 成立。

进一步地， $(mp)\cdot(mr) = (x_{1}^{2}+x_{2}^{2}+x_{3}^{2}+x_{4}^{2})\cdot(y_{1}^{2}+y_{2}^{2}+y_{3}^{2}+y_{4}^{2}) = z_{1}^{2}+z_{2}^{2}+z_{3}^{2}+z_{4}^{2}$ ，这里的 $z_{i}$ 正如恒定式里面所定义的。由于 $y_{i}\equiv x_{i}(\mod m)$ ，并且 $\sum_{i=1}^{4}x_{i}^{2}\equiv 0(\mod m)$ 。因此， $z_{i}\equiv 0(\mod m)$ 对于 $i\in\{1,2,3,4\}$ 都成立。所以， $z_{i}=w_{i}m$ ， $w_{i}\in\mathbb{Z}$ 对于 $i\in\{1,2,3,4\}$ 都成立。通过 $(mp)\cdot(mr) = \sum_{i=1}^{4}z_{i}^{2}$ 可以得到 $pr=\sum_{i=1}^{4}w_{i}^{2}$ 成立。但是， $1\leq r<m$ 这与 $m$ 的最小性假设矛盾了。

因此， $m=1$ ，Claim 证明完毕。

于是，对于所有的奇素数，都可以表示为四个整数的平方之和。根据之前的分析，可以得到对于所有的正整数，都可以表示为四个整数的平方之和。Lagrange 定理证明完毕。

参考文献

GTM 164, Additive Number Theory, Melvyn B.Nathanson, 1996.

Mathematics

传染病的数学模型

January 25, 2020 zr9558 Leave a comment

近期，国内的疫情闹得沸沸扬扬，很多省市自治区都出现了流感的患者。回想起之前在学校的时候曾经研究过微分方程和动力系统，于是整理一下相关的数学模型，分享给各位读者。笔者并不是研究这个领域的专家，并且这篇文章只是从微分方程角度出发，分析方程的性质，不一定适用于真实环境，而且真实环境比这个也复杂得多。

关于传染病的数学模型，在许多年前数学界早已做过研究，根据传染病的传播速度不同，空间范围各异，传播途径多样，动力学机理等各种因素，对传染病模型按照传染病的类型划分为 SI，SIR，SIRS，SEIR 模型。如果是按照连续时间来划分，那么这些模型基本上可以划分为常微分方程（Ordinary Differential Equation），偏微分方程（Partial Differential Equation）等多种方程模型；如果是基于离散的时间来划分，那么就是所谓的差分方程（Difference Equation）。

在本文中，将会主要介绍常微分方程中的一些传染病数学模型。在介绍方程之前，首先要介绍一些常用的符号：在时间戳 $t$ 上，可以定义以下几种人群：

易感者（susceptible）：用符号 $S(t)$ 来表示；
感染者（infective）：用符号 $I(t)$ 来表示；
康复者（Recoverd）：用符号 $R(t)$ 来表示；

其次，在时间戳 $t$ 上，总人口是 $N(t)=S(t)+I(t)+R(t)$ 。如果暂时不考虑人口增加和死亡的情况，那么 $N(t)\equiv N$ 是一个恒定的常数值。

除此之外，

$r$ 表示在单位时间内感染者接触到的易感者人数；
传染率： $\beta$ 表示感染者接触到易感者之后，易感者得病的概率；
康复率： $\gamma$ 表示感染者康复的概率，有可能变成易感者（可再感染），也有可能变成康复者（不再感染）。

在进行下面的分析之前，先讲一个常微分方程的解。

Claim. 假设 $x=x(t)$ 是关于 $t$ 的一个方程，且满足 $\frac{dx}{dt} + a_{1}x + a_{2}x^{2}=0$ 和 $x(0)=x_{0}$ ，那么它的解是： $x(t) = \frac{e^{-a_{1}t}}{\frac{1}{x_{0}}-\frac{a_{2}}{a_{1}}(e^{-a_{1}t}-1)}$ .

Proof. 证明如下：

通过 $\frac{dx}{dt}+a_{1}x+a_{2}x^{2}=0$ 可以得到 $-\frac{d}{dt}\bigg(\frac{1}{x}\bigg) + a_{1}\bigg(\frac{1}{x}\bigg)+a_{2}=0$ ；令 $y = 1/x$ ，得到 $\frac{dy}{dt}-a_{1}y=a_{2}$ 。所以， $\frac{d}{dt}(e^{-a_{1}t}y) = a_{2}e^{-a_{1}t}$ ，两边积分可以得到 $e^{-a_{1}t}y-y_{0}=\bigg(-\frac{a_{2}}{a_{1}}\bigg)(e^{-a_{1}t}-1)$ ，其中 $y_{0}=1/x_{0}$ 。求解之后得到： $x(t) = e^{-a_{1}t}/\bigg(\frac{1}{x_{0}}-\frac{a_{2}}{a_{1}}(e^{-a_{1}t}-1)\bigg)$ 。

SI 模型（Susceptible-Infective Model）

在 SI 模型里面，只考虑了易感者和感染者，并且感染者不能够恢复，此类病症有 HIV 等；

其微分方程就是：

$\begin{cases}\frac{dS}{dt} = -\frac{r\beta I}{N} S \\ \frac{dI}{dt}=\frac{r\beta I}{N}S \end{cases}$

初始条件就是 $S(0)=S_{0}$ ， $I(0) = I_{0}$ ，并且 $S(t)+I(t)=N$ 对于所有的 $t\geq 0$ 都成立。

于是，把 $S = N - I$ 代入第二个微分方程可以得到： $\frac{dI}{dt} - r\beta I + \frac{r\beta}{N}I^{2}=0$ 。因此根据前面所提到的常微分方程的解可以得到：

$I(t) = \frac{NI_{0}}{I_{0}+(N-I_{0})e^{-r\beta t}}$ .

这个就是所谓的逻辑回归函数，而在机器学习领域，最简单的逻辑回归函数就是 $\sigma(x) = 1/(1+e^{-x})$ 这个定义。而 $I(t)$ 只是做了一些坐标轴的平移和压缩而已。由于 $\lim_{t\rightarrow +\infty}e^{-t}=0$ ，所以， $\lim_{t\rightarrow +\infty}I(t) = N$ ，从而 $\lim_{t\rightarrow +\infty}S(t) = 0$ 。

通过数值模拟可以进一步知道：

简单来看，在 SI 模型的假设下，全部人群到最后都会被感染。

SIS 模型（Susceptible-Infectious-Susceptible Model）

除了 HIV 这种比较严重的病之外，还有很多小病是可以恢复并且反复感染的，例如日常的感冒，发烧等。在这种情况下，感染者就有一定的几率重新转化成易感者。如下图所示：

其微分方程就是：

$\begin{cases} \frac{dS}{dt} = -r \beta S\frac{I}{N} + \gamma I \\ \frac{dI}{dt}=r\beta S \frac{I}{N} - \gamma I \end{cases}$ ，其初始条件就是 $S(0)=S_{0}$ ， $I(0)=I_{0}$ .

使用同样的方法，把 $S=N-I$ 代入第二个微分方程可以得到： $\frac{dI}{dt} - (r\beta - \gamma)I + \frac{r\beta}{N}I^{2}=0$ . 通过之前的 Claim 可以得到解为：

$I(t) = \frac{N(r\beta-\gamma)}{r\beta}/\bigg(\bigg(\frac{N(r\beta-\gamma)}{I_{0}r\beta}-1\bigg)e^{-(r\beta-\gamma)t}+1\bigg)$ .

从而可以得到 $\lim_{t\rightarrow +\infty} I(t) = N(r\beta - \gamma)/(r\beta)$ 且 $\lim_{t\rightarrow +\infty} S(t) = (N\gamma)/(r\beta)$ . 这个方程同样也是逻辑回归方程，只是它的渐近线与之前的 SI 模型有所不同。

SIR 模型（Susceptible-Infectious-Recovered Model）

有的时候，感染者在康复了之后，就有了抗体，于是后续就不再会获得此类病症，这种时候，考虑 SIS 模型就不合适了，需要考虑 SIR 模型。此类病症有麻疹，腮腺炎，风疹等。

其微分方程是： $\begin{cases} \frac{dS}{dt}=-r\beta S \frac{I}{N} \\ \frac{dI}{dt}=r\beta S\frac{I}{N} - \gamma I \\ \frac{dR}{dt}=\gamma I\end{cases}$ 。其初始条件是 $S(0)=S_{0}, I(0)=I_{0}, R(0)=R_{0}$ ，并且 $S(t), I(t), R(t)\geq 0$ 和 $S(t) +I(t)+R(t)=N$ 对于所有的 $t\geq 0$ 都成立。

对于这类方程，就不能够得到其解析解了，只能够从它的动力系统开始进行分析，得到解的信息。根据第一个微分方程可以得到： $\frac{dS}{dt} = -r\beta S\frac{I}{N}<0$ ，于是 $S(t)$ 是一个严格递减函数。同时， $0\leq S(t)\leq N$ 对于所有的 $t\geq 0$ 都成立，于是存在 $S_{\infty}\in[0,\infty]$ 使得 $\lim_{t\rightarrow \infty}S(t)=S_{\infty}$ .

通过第一个微分方程和第二个微分方程可以得到： $\frac{d(S+I)}{dt} = -\gamma I$ ，因此对它两边积分得到 $\int_{0}^{T} \frac{d(S+I)}{dt} = -\gamma \int_{0}^{T}I(t)dt$ . 左侧等于 $S(T) + I(T) - S(0) - I(0)$ ，上界是 $4N$ ，因此令 $T\rightarrow \infty$ 可以得到 $\int_{0}^{\infty}I(t) dt\leq 4N/\gamma$ . 而 $I(t)\geq 0$ 且是连续可微函数，因此 $\lim_{t\rightarrow \infty}I(t) = 0$ 。这意味着所有的感染人群都将康复。

由于 $S(t)$ 是严格单调递减函数，因此从第二个微分方程可以得到：当 $S(t) = N\gamma/(r\beta)$ 时，感染人数 $I(t)$ 达到最大值。

其余模型

在以上的 SI，SIS，SIR 模型中，还可以把死亡因素考虑进去。除此之外，还有 SIRS 模型，SEIR 模型等，在这里就不再做赘述。有兴趣的读者可以参阅相关的参考书籍。

参考文献

Introduction to SEIR Models, Nakul Chitnis, Workshop on Mathematical Models of Climate Variability, Environmental Change and Infectious Diseases, Trieste, Italy, 2017

数据挖掘与机器学习

符号计算中的深度学习方法

January 11, 2020 zr9558 Leave a comment

符号计算

符号计算一直是计算数学的重要领域之一。在开源领域，Python 的 SymPy 就可以支持符号计算。在商业化领域，Maple，Matlab，Mathematica 都能够进行符号计算。它们不仅能够做简单的实数和复数加减乘除，还能够支持数学分析，线性代数，甚至各种各样的大学数学课程。

随着人工智能的进一步发展，深度学习不仅在图像识别，自然语言处理方向上发挥着自身的价值，还在各种各样的领域展示着自己的实用性。在 2019 年底，facebook 两位研究员在 arxiv 上挂出了一篇文章《Deep Learning for Symbolic Mathematics》，在符号计算方向上引入了深度学习的工具。

要想了解符号运算，就要先知道在计算机中，是怎么对数学公式进行表示的。较为常见的表达式形如：

$2 + 3 * (5 + 2)$
$3x^{2}+\cos(2x)+1$
$\frac{\partial^{2}\psi}{\partial x^{2}} - \frac{1}{v^{2}}\frac{\partial^{2}\psi}{\partial t^{2}}$

在这里，数学表达式通常都会被表示成树的结构，其中树的内部节点是由算子（operator），函数（function）组成的，叶子节点由数字，变量，函数等组成。例如：

图 1 的三幅图分别对应着上面的三个数学表达式。

在 Python 的 SymPy 工具中，同样可以对数学公式进行展示。其表示方法就是用 sympy.srepr

>>> import sympy
>>> x, y = sympy.symbols("x y")
>>> expr = sympy.sin(x+y) + x**2 + 1/y - 10
>>> sympy.srepr(expr)
"Add(Pow(Symbol('x'), Integer(2)), sin(Add(Symbol('x'), Symbol('y'))), Integer(-10), Pow(Symbol('y'), Integer(-1)))"
>>> expr = sympy.sin(x*y)/2 - x**2 + 1/y
>>> sympy.srepr(expr)
"Add(Mul(Integer(-1), Pow(Symbol('x'), Integer(2))), Mul(Rational(1, 2), sin(Mul(Symbol('x'), Symbol('y')))), Pow(Symbol('y'), Integer(-1)))"

SymPy 的 srepr 函数的输出用树状结构来表示就是形如图 2 这种格式。叶子节点要么是 x，y 这种变量，要么是 -1 和 2 这种整数。对于一元函数而言，例如 sin 函数，就是对应唯一的一个叶子。对于二元函数而言，例如 pow，mul，add，则是对应两个叶子节点。

论文方案

在 Deep Learning for Symbolic Mathematics 这篇论文中，作者们的大致思路是分成以下几步的：

生成数据；
训练模型；
预测结果；

第一步生成数据是为了让深度学习模型是大量的已知样本来做训练；第二步训练模型是用第一步的数据把端到端的深度学习模型进行训练；第三步预测结果是给一个函数或者一个微分方程，使用已经训练好的模型来预测结果，对预测出来的多个结果进行排序，选择最有可能的那个结果作为符号计算的值。

众所周知，深度学习的训练是依赖大量的样本数据的，那么要想用深度学习来解决符号计算的问题，就要解决样本少的问题。在这篇论文中，作者们把精力投入了三个领域，分别是：

函数的积分；
一阶常微分方程；
二阶常微分方程。

在生成数据之前，作者们对数据的范围进行了必要的限制：

数学表达式最多拥有 15 个内部节点；
$L = 11$ 表示叶子节点的值只有 11 个，分别是变量 $x$ 和 $\{-5,-4,-3,-2,-1,1,2,3,4,5\}$ ；
$p_{1} = 15$ 表示一元计算只有 15 个，分别是 $\exp, \log, \sqrt, \sin, \cos, \tan, \arcsin, \arccos, \arctan$ , $sinh, cosh, tanh, arcsinh, arccosh, arctanh$ 。
$p_{2} = 4$ 表示二元计算只有四个，分别是 +, -, *, /；

意思就是在这个有限的范围内去生成深度学习所需要的数据集。

积分数据的生成

在微积分里面，积分指的是求导的逆运算，那么形如 $(f',f)$ 的表达式就可以作为深度学习的积分训练数据。生成积分的话其实有多种方法：

第一种方法：前向生成（Forward Generation，简写为 FWD）。主要思路就是在以上的数据范围内随机生成各种各样的方程 $f$ ，然后使用 SymPy 或者 Mathematica 等工具来计算函数 $f$ 的积分 $F$ ，那么 $(f,F)$ 就可以作为一个训练集。当然，有的时候函数 $f$ 的积分是无法计算出来的，那么这种计算表达式就需要进行放弃，就不能放入训练集。

第二种方法：反向生成（Backward Generation，简写为 BWD）。由于积分是求导的逆运算，可以在以上的数据范围内随机生成各种各样的方程 $f$ ，然后计算它们的导数 $f'$ ，于是 $(f',f)$ 就可以放入积分数据的训练集。

第三种方法：分部积分（Backward generation with integration by parts，简写为 IBP）。根据分部积分的公式，如果 $F'=f, G'=g$ ，那么 $\int F(x)g(x)dx=F(x)G(x) - \int f(x)G(x)dx$ 。对于两个随机生成的函数 $F, G$ ，可以计算出它们的导数 $f, g$ 。如果 $fG$ 在训练集合里面，那么就把 $Fg$ 的积分计算出来放入训练集合；反之，如果 $Fg$ 在训练集合里面，那么就把 $fG$ 的积分计算出来放入训练集合。如果 $fG$ 和 $Fg$ 都没有在训练集合，并且都无法积分出来，那么就放弃该条数据。

三种方法的比较可以参见表 1，从表 1 可以看出通过分部积分（Integration by parts）所获得的样本表达式相对于前向生成和反向生成都比较长。

一阶常微分方程的生成

一阶常微分方程只有该函数的一阶导数，因此在构造函数的时候，作者们用了一个技巧。在随机生成表达式的时候，叶子节点的元素只从 $\{x, -5,-4,-3,-2,-1,1,2,3,4,5\}$ 选择一个，于是随机把其中的一个整数换成变量 $c$ 。例如：在 $x\log(2/x)$ 中就把 2 换成 c，于是得到了一个二元函数 $f(x,c) = x\log(c/x)$ 。那么就执行以下步骤：

生成二元函数 $f(x,c) = x\log(c/x)$ ；
求解 $c$ ，得到 $c = xe^{f(x)/x}$ ；
对 $x$ 进行求导得到 $0 = e^{f(x)/x}(1+f'(x)-f(x)/x) = 0$ ；
简化后得到 $x+xf'(x) - f(x) =0$ ，也就是 $x+xy'-y=0$ 。

此时，一阶常微分方程的训练数据就是 $(x+xy'-y=0, x\log(c/x))$ 。

二阶常微分方程不仅可能由该函数的一阶导数，还必须有二阶导数，那么此时就要求导两次：

生成三元函数 $f(x,c_{1},c_{2}) = c_{1}e^{x}+c_{2}e^{-x}$ ；
求解 $c_{2}$ 得到 $c_{2} = f(x,c_{1},c_{2}) e^{x}- c_{1}e^{2x}$ ；
对 $x$ 进行求导得到 $0 = e^{x}(\partial f(x,c_{1},c_{2})/\partial x + f(x,c_{1},c_{2})) - 2c_{1}e^{2x} = 0$ ；
求解 $c_{1}$ 得到 $c_{1} = e^{-x}(\partial f(x,c_{1},c_{2})/\partial x+f(x))/2$ ；
对 $x$ 进行求导得到 $0 = e^{-x}(\partial^{2} f(x,c_{1},c_{2})/\partial x^{2} - f(x,c_{1},c_{2}))=0$ ；
简化后得到 $\partial^{2} f(x,c_{1},c_{2})/\partial x^{2} - f(x,c_{1},c_{2})=0$ ，也就是 $y''-y=0$ 。

那么此时的二阶微分方程的训练数据就是 $(y''-y=0, c_{1}e^{x}+c_{2}e^{-x})$ 。

需要注意的事情就是，在生成数据的时候，一旦无法求出 $c_{1}, c_{2}$ 的表达式，该数据就需要放弃，重新生成新的数据。

数据处理

数学表达式的简化（expression simplification）：例如 $x+1+1$ 可以简化成 $x+2$ ， $\sin^{2}(x)+\cos^{2}(x)$ 可以简化成 1。
参数的简化（coefficient simplification）：例如 $\log(x^{2}) + c\log(x)$ 可以简化成 $c\log(x)$ 。
无效表达式的过滤（invalid expression filter）：例如 $\sqrt{2}, \log(0)$ 等。

树状结构的表达式，是使用前缀表达式来写成一个句子的。例如 2+3 就写成 + 2 3，2 + x 就写成 + 2 x。

模型训练

在这里，作者们用了 Transformer 模型，8 attention heads，6 layers，512 个维度，并且发现用更复杂的模型并没有提升其效果。在预测的时候，使用不同的 Beam Size，其准确率是不一样的，在 Beam Size = 50 的时候，效果比较好。参见表 2。

在与其他数学软件的比较中，作者限制了 Mathematica 的运行时间为 30s。

并且举出了几个现有模型能够计算出来，Mathematica 无法计算的例子。

在求解的等价性方面，可以根据 Score 逆序排列，然后这个例子的 Top10 都是等价的。

在使用 SymPy 的过程中，可以获得各种各样的积分表达式如下：

结论

符号计算在 1960 年代末就已经在研究了，有诸多的符号计算软件，例如 Matlab，Mathematica，Maple，PARI，SAGE 等。在这篇论文中，作者们使用标准的 seq2seq 模型来对生成的数据进行训练，然后获得积分，一阶常微分方程，二阶常微分方程的解。在传统符号计算的基础上，开拓了一种新的思路。

Mathematics

用 Python 来研究数学 — SymPy 符号工具包介绍

December 12, 2019 zr9558 Leave a comment

SymPy 的简单介绍

SymPy 是一个符号计算的 Python 库，完全由 Python 写成，为许多数值分析，符号计算提供了重要的工具。SymPy 的第一个版本于 2007 年开源，并且经历了十几个版本的迭代，在 2019 年已经基于修正的 BSD 许可证开源了 1.4 版本。SymPy 的开源地址和官方网站分别是：

GitHub 链接：https://github.com/sympy/sympy
SymPy 官方网站：https://www.sympy.org/en/index.html

从 SymPy 的 1.4 版本文档中，可以看出，SymPy 可以支持很多初等数学，高等数学，甚至研究生数学的符号计算。在初等数学和高等数学中，SymPy 可以支持的内容包括但不限于：

基础计算（Basic Operations）；
公式简化（Simplification）；
微积分（Calculus）；
解方程（Solver）；
矩阵（Matrices）；
几何（geometry）；
级数（Series）；

在更多的数学领域中，SymPy 可以支持的内容包括但不限于：

范畴论（Category Theory）；
微分几何（Differential Geometry）；
常微分方程（ODE）；
偏微分方程（PDE）；
傅立叶变换（Fourier Transform）；
集合论（Set Theory）；
逻辑计算（Logic Theory）。

SymPy 的工具库介绍

SymPy 的基础计算

在数学中，基础的计算包括实数和复数的加减乘除，那么就需要在程序中描述出实数与复数。著名的欧拉公式

$e^{i\pi}+1 = 0$ 。

正好用到了数学中最常见的五个实数。在 SymPy 里面， $e, i, \pi, \infty$ 是用以下符号来表示的：其中 sympy.exp() 表示以 $e$ 为底的函数。

sympy.exp(1), sympy.I, sympy.pi, sympy.oo

而想要计算欧拉公式的话，只需要输入下面的公式即可：

>>> sympy.exp(sympy.I * sympy.pi) + 1
0

如果需要看 $e, \pi$ 的小数值，可以使用 evalf() 函数，其中 evalf() 函数里面的值表示有效数字的位数。例如下面就是精确到 10 位有效数字。当然，也可以不输入。

>>> sympy.E.evalf(10)
2.718281828
>>> sympy.E.evalf()
2.71828182845905
>>> sympy.pi.evalf(10)
3.141592654
>>> sympy.pi.evalf()
3.14159265358979

除此之外，如果需要查看某个实数的有效数字，也是类似操作的：

>>> expr = sympy.sqrt(8)
>>> expr.evalf()
2.82842712474619

而对于实数的加减乘除，则可以如下操作：

>>> x, y= sympy.symbols("x y")
>>> x + y
x + y
>>> x - y
x - y
>>> x * y
x*y
>>> x / y
x/y

而对于复数的加减乘除，则是类似的操作，令两个复数分别是 $z_{1} = x_{1} + i y_{1}$ 和 $z_{2} = x_{2} + i y_{2}$ 。

>>> x1, y1, x2, y2 = sympy.symbols("x1 y1 x2 y2")
>>> z1 = x1 + y1 * sympy.I
x1 + I*y1
>>>  z2 = x2 + y2 * sympy.I
x2 + I*y2
>>> z1 + z2
x1 + x2 + I*y1 + I*y2
>>> z1 - z2
x1 - x2 + I*y1 - I*y2
>>> z1 * z2
(x1 + I*y1)*(x2 + I*y2)
>>> z1 / z2
(x1 + I*y1)/(x2 + I*y2)

对于多项式而言，有的时候我们希望将其展开，有的时候则需要将其合并，最终将其简化成最简单的形式。

>>> sympy.expand((x+1)**2)
x**2 + 2*x + 1
>>> sympy.expand((x+1)**5)
x**5 + 5*x**4 + 10*x**3 + 10*x**2 + 5*x + 1
>>> sympy.factor(x**3+1)
(x + 1)*(x**2 - x + 1)
>>> sympy.factor(x**2+3*x+2)
(x + 1)*(x + 2)
>>> sympy.simplify(x**2 + x + 1 - x)
x**2 + 1
>>> sympy.simplify(sympy.sin(x)**2 + sympy.cos(x)**2)
1

在多变量的场景下，SymPy 也可以对其中的某个变量合并同类项，同时还可以计算某个变量的某个次数所对应的系数是多少，例如：

>>> expr = x*y + x - 3 + 2*x**2 - x**2 + x**3 + y**2 + x**2*y**2
>>> sympy.collect(expr,x)
x**3 + x**2*(y**2 + 1) + x*(y + 1) + y**2 - 3
>>> sympy.collect(expr,y)
x**3 + x**2 + x*y + x + y**2*(x**2 + 1) - 3
>>> expr.coeff(x, 2)
y**2 + 1
>>> expr.coeff(y, 1)
x

有理函数形如 $f(x) = p(x)/q(x)$ ，其中 $p(x)$ 和 $q(x)$ 都是多项式。一般情况下，我们希望对有理函数进行简化，合并或者分解的数学计算。

在需要合并的情形下，如果想把有理函数处理成标准格式 $p(x)/q(x)$ 并且去除公因子，那么可以使用 cancel 函数。另一个类似的就是 together 函数，但是不同之处在于 cancel 会消除公因子，together 不会消除公因子。例如：

$expr = \frac{x^{2}+3x+2}{x^{2}+x}$ ，

>>> expr = (x**2 + 3*x + 2)/(x**2 + x)
>>> sympy.cancel(expr)
(x + 2)/x
>>> sympy.together(expr)
(x**2 + 3*x + 2)/(x*(x + 1))

除了合并和消除公因子之外，有的时候还希望对分子和分母进行因式分解，例如：

expr = (x**2 + 3*x + 2)/(x**2 + x)
>>> sympy.factor(expr)
(x + 2)/x
>>> expr = (x**3 + 3*x**2 + 2*x)/(x**5+x)
>>> sympy.factor(expr)
(x + 1)*(x + 2)/(x**4 + 1)
>>> expr = x**2 + (2*x+1)/(x**3+1)
>>> sympy.factor(expr)
(x**5 + x**2 + 2*x + 1)/((x + 1)*(x**2 - x + 1))

合并的反面就是部分分式展开（Partial Fraction Decomposition），它是把有理函数分解成多个次数较低的有理函数和的形式。这里需要用 apart 函数：

>>> expr = (x**4 + 3*x**2 + 2*x)/(x**2+x)
>>> sympy.apart(expr)
x**2 - x + 4 - 2/(x + 1)
>>> expr = (x**5 + 1)/(x**3+1)
>>> sympy.apart(expr)
x**2 - (x - 1)/(x**2 - x + 1)

在 SymPy 里面，同样支持各种各样的三角函数，例如：三角函数的简化函数 trigsimp，三角函数的展开 expand_trig，

>>> expr = sympy.sin(x)**2 + sympy.cos(x)**2
>>> sympy.trigsimp(expr)
1
>>> sympy.expand_trig(sympy.sin(x+y))
sin(x)*cos(y) + sin(y)*cos(x)
>>> sympy.expand_trig(sympy.cos(x+y))
-sin(x)*sin(y) + cos(x)*cos(y)
>>> sympy.trigsimp(sympy.sin(x)*sympy.cos(y) + sympy.sin(y)*sympy.cos(x))
sin(x + y)
>>> sympy.trigsimp(-sympy.sin(x)*sympy.sin(y) + sympy.cos(x)*sympy.cos(y))
cos(x + y)

同样的，在乘幂上面，同样有简化函数 powsimp，效果与之前提到的 simplify 一样。除此之外，还可以根据底数来做合并，即分别使用 expand_power_exp 函数与 expand_power_base 函数。

>>> sympy.powsimp(x**z*y**z*x**z)
x**(2*z)*y**z
>>> sympy.simplify(x**z*y**z*x**z)
x**(2*z)*y**z
>>> sympy.expand_power_exp(x**(y + z))
x**y*x**z
>>> sympy.expand_power_base(x**(y + z))
x**(y + z)

作为指数的反函数对数，sympy 也是有着类似的展开合并函数，expand_log，logcombine 承担了这样的角色。

$\ln(xy) = \ln(x) + \ln(y)$ ，

$\ln(x/y) = \ln(x) - \ln(y)$ 。

>>> sympy.expand_log(sympy.log(x*y), force=True)
log(x) + log(y)
>>> sympy.expand_log(sympy.log(x/y), force=True)
log(x) - log(y)

SymPy 的微积分工具

下面，我们会从一个最基本的函数 $f(x) = 1/x$ 出发，来介绍 SymPy 的各种函数使用方法。如果想进行变量替换，例如把 $x$ 变成 $y$ ，那么可以使用 substitution 方法。除此之外，有的时候也希望能够得到函数 $f$ 在某个点的取值，例如 $f(1) = 1/1 = 1$ ，那么可以把参数换成 1 即可得到函数的取值。例如，

>>> import sympy
>>> x = sympy.Symbol("x")
>>> f = 1 / x
1/x
>>> y = sympy.Symbol("y")
>>> f = f.subs(x,y)
1/y
>>> f = f.subs(y,1)
1

在微积分里面，最常见的概念就是极限，SymPy 里面的极限函数就是 limit。使用方法如下：

>>> f = 1/x
>>> sympy.limit(f,x,0)
oo
>>> sympy.limit(f,x,2)
1/2
>>> sympy.limit(f,x,sympy.oo)
0
>>> g = x * sympy.log(x)
>>> sympy.limit(g,x,0)
0

对于函数 $f(x) = 1/x$ 而言，它的导数计算函数是 diff， $n$ 阶导数也可以用这个函数算。

>>> f = 1/x
>>> sympy.diff(f,x)
-1/x**2
>>> sympy.diff(f,x,2)
2/x**3
>>> sympy.diff(f,x,3)
-6/x**4
>>> sympy.diff(f,x,4)
24/x**5

提到 $n$ 阶导数，就必须要提一下 Taylor Series 了。对于常见函数的 Taylor Series，SymPy 也是有非常简便的方法，那就是 series 函数。其参数包括 expr, x, x0, n, dir，分别对应着表达式，函数的自变量，Taylor Series 的中心点，n 表示阶数，dir 表示方向，包括”+-“，”-“，”+”，分别表示 $x\rightarrow x0$ , $x\rightarrow x0^{-}$ , $x\rightarrow x0^{+}$ 。

sympy.series.series.series(expr, x=None, x0=0, n=6, dir='+') >>> g = sympy.cos(x) >>> sympy.series(g, x) 1 - x**2/2 + x**4/24 + O(x**6) >>> sympy.series(g, x, x0=1, n=10) cos(1) - (x - 1)*sin(1) - (x - 1)**2*cos(1)/2 + (x - 1)**3*sin(1)/6 + (x - 1)**4*cos(1)/24 - (x - 1)**5*sin(1)/120 - (x - 1)**6*cos(1)/720 + (x - 1)**7*sin(1)/5040 + (x - 1)**8*cos(1)/40320 - (x - 1)**9*sin(1)/362880 + O((x - 1)**10, (x, 1))

积分的计算函数是 integrate，包括定积分与不定积分：

$\int\frac{1}{x}dx = \ln(x)+C$ ，

$\int_{1}^{2}\frac{1}{x}dx = \ln(2)$ 。

>>> f = 1/x
>>> sympy.integrate(f,x)
log(x)
>>> sympy.integrate(f, (x,1,2))
log(2)

对于广义积分而言，就需要用到 $\infty$ 这个概念了，但是在 SymPy 里面的写法还是一样的。

$\int_{-\infty}^{0}e^{-x^{2}}dx=\frac{\sqrt{\pi}}{2}$ ，

$\int_{0}^{+\infty}e^{-x}dx = 1$ ，

$\int_{-\infty}^{+\infty}\int_{-\infty}^{+\infty}e^{-x^{2}-y^{2}}dxdy = \pi$ 。

>>> g = sympy.exp(-x**2)
>>> sympy.integrate(g, (x,-sympy.oo,0))
sqrt(pi)/2
>>> g = sympy.exp(-x)
>>> sympy.integrate(g, (x, 0, sympy.oo))
1
>>> h = sympy.exp(-x**2 - y**2)
>>> sympy.integrate(h, (x,-sympy.oo, sympy.oo), (y, -sympy.oo, sympy.oo))
pi

SymPy 的方程工具

在初等数学中，通常都存在一元一次方程，一元二次方程等，并且在不同的域上有着不同的解。SymPy 里面的相应函数就是 solveset，根据定义域的不同，可以获得完全不同的解。

$\{x\in\mathbb{R}: x^{3}-1=0\}$ ，

$\{x\in\mathbb{C}:x^{3}-1=0\}$ ，

$\{x\in\mathbb{R}:e^{x}-x=0\}$ ，

$\{x\in\mathbb{R}:e^{x}-1=0\}$ ，

$\{x\in\mathbb{C}:e^{x}-1=0\}$ 。

>>> sympy.solveset(sympy.Eq(x**3,1), x, domain=sympy.S.Reals)
{1}
>>> sympy.solveset(sympy.Eq(x**3,1), x, domain=sympy.S.Complexes)
{1, -1/2 - sqrt(3)*I/2, -1/2 + sqrt(3)*I/2}
>>> sympy.solveset(sympy.Eq(x**3 - 1,0), x, domain=sympy.S.Reals)
{1}
>>> sympy.solveset(sympy.Eq(x**3 - 1,0), x, domain=sympy.S.Complexes)
{1, -1/2 - sqrt(3)*I/2, -1/2 + sqrt(3)*I/2}
>>> sympy.solveset(sympy.exp(x),x)
EmptySet()
>>> sympy.solveset(sympy.exp(x)-1,x,domain=sympy.S.Reals)
{0}
>>> sympy.solveset(sympy.exp(x)-1,x,domain=sympy.S.Complexes)
ImageSet(Lambda(_n, 2*_n*I*pi), Integers)

在这里，Lambda 表示的是数学公式，第一个是自变量，第二个是函数，最后是自变量的定义域。

在线性代数中，最常见的还是多元一次方程组，那么解法是一样的：

$\begin{cases}x+y-10=0 \\ x-y-2=0\end{cases}$

>>> sympy.solve([x+y-10, x-y-2], [x,y])
{x: 6, y: 4}

对于三角函数，也是类似的写法：

$\begin{cases} \sin(x-y)=0 \\ \cos(x+y)=0 \end{cases}$ ，

>>> sympy.solve([sympy.sin(x-y), sympy.cos(x+y)], [x,y])
[(-pi/4, 3*pi/4), (pi/4, pi/4), (3*pi/4, 3*pi/4), (5*pi/4, pi/4)]

SymPy 的矩阵工具

在矩阵论中，最常见的就是单位矩阵了，而单位矩阵只与一个参数有关，那就是矩阵的大小。下面就是 3*3，3*2，2*3 大小的矩阵。

>>> sympy.eye(3)
Matrix([
[1, 0, 0],
[0, 1, 0],
[0, 0, 1]])
>>> sympy.eye(3,2)
Matrix([
[1, 0],
[0, 1],
[0, 0]])
>>> sympy.eye(2,3)
Matrix([
[1, 0, 0],
[0, 1, 0]])

另外还有全零和全一矩阵，意思就是矩阵中的所有值全部是零和一。

>>> sympy.ones(2,3)
Matrix([
[1, 1, 1],
[1, 1, 1]])
>>> sympy.zeros(3,2)
Matrix([
[0, 0],
[0, 0],
[0, 0]])

而对角矩阵也可以使用 diag 轻松获得：

>>> sympy.diag(1,1,2)
Matrix([
[1, 0, 0],
[0, 1, 0],
[0, 0, 2]])

而矩阵的加法，减法，乘法，逆运算，转置，行列式，SymPy 都是可以支持的：

$A = \left(\begin{array}{cc} 1 & 1 \\ 0 & 2 \end{array}\right)$ ，

$B = \left(\begin{array}{cc} 1 & 0 \\ 1 & 1 \end{array}\right)$ ，

>>> A = sympy.Matrix([[1,1],[0,2]])
>>> B = sympy.Matrix([[1,0],[1,1]])
>>> A
Matrix([
[1, 1],
[0, 2]])
>>> B
Matrix([
[1, 0],
[1, 1]])
>>> A + B
Matrix([
[2, 1],
[1, 3]])
>>> A - B
Matrix([
[ 0, 1],
[-1, 1]])
>>> A * B
Matrix([
[2, 1],
[2, 2]])
>>> A * B**-1
Matrix([
[ 0, 1],
[-2, 2]])
>>> B**-1
Matrix([
[ 1, 0],
[-1, 1]])
>>> A.T
Matrix([
[1, 0],
[1, 2]])
>>> A.det()
2

在某些情况下，需要对矩阵进行加上一行或者加上一列的操作，在这里就可以使用 row_insert 或者 col_insert 函数：第一个参数表示插入的位置，第二个参数就是相应的行向量或者列向量。而在删除上就很简单了，直接使用 row_del 或者 col_del 即可。

>>> A
Matrix([
[1, 1],
[0, 2]])
>>> A.row_insert(1, sympy.Matrix([[10,10]]))
Matrix([
[ 1, 1],
[10, 10],
[ 0, 2]])
>>> A.col_insert(0, sympy.Matrix([3,3]))
Matrix([
[3, 1, 1],
[3, 0, 2]])
>>> A.row_del(0)
>>> A
Matrix([[0, 2]])
>>> A.col_del(1)
>>> A
Matrix([[0]])

在对角化方面，同样可以使用 eigenvals()，eigenvecs(), diagonalize() 函数：

>>> A
Matrix([
[1, 1],
[0, 2]])
>>> A.eigenvals()
{2: 1, 1: 1}
>>> A.eigenvects()
[(1, 1, [Matrix([
[1],
[0]])]), (2, 1, [Matrix([
[1],
[1]])])]
>>> A.diagonalize()
(Matrix([
[1, 1],
[0, 1]]), Matrix([
[1, 0],
[0, 2]]))

在 eigenvals() 返回的结果中，第一个表示特征值，第二个表示该特征值的重数。在特征向量 eigenvecs() 中，第一个表示特征值，第二个表示特征值的重数，第三个表示特征向量。在对角化 diagonalize() 中，第一个矩阵表示 $P$ ，第二个矩阵表示 $D$ ， $A = P*D*P^{-1}$ 。

在矩阵中，最常见的还是多元一次方程的解。如果要求 $Ax =b$ 的解，可以有以下方案：

>>> A = sympy.Matrix([[1,1],[0,2]])
>>> A
Matrix([
[1, 1],
[0, 2]])
>>> b = sympy.Matrix([3,5])
>>> b
Matrix([
[3],
[5]])
>>> A**-1*b
Matrix([
[1/2],
[5/2]])
>>> sympy.linsolve((A,b))
{(1/2, 5/2)}
>>> sympy.linsolve((A,b),[x,y])
{(1/2, 5/2)}

SymPy 的集合论工具

集合论可以说是数学的基础，在任何数学的方向上都能够看到集合论的身影。在 SymPy 里面，有一个类叫做 sympy.sets.sets.set。在集合论里面，常见的就是边界，补集，包含，并集，交集等常见的操作。但是感觉 SymPy 中的集合论操作主要集中在实数域或者复数域上。

对于闭区间 $I=[0,1]$ 和开区间 $J = (0,1)$ 而言，在 SymPy 中使用以下方法来表示：

I = sympy.Interval(0,1)
J = sympy.Interval.open(0,1)
K = sympy.Interval(0.5,2)

其开始和结束的点可以分别使用 start 和 end 来表示：

>>> I.start
0
>>> I.end
1

其长度用 measure 来表示：

>>> I.measure
1

其闭包用 closure 来表示：

>>> I.closure
Interval(0, 1)

其内点用 interior 来表示：

>>> I.interior
Interval.open(0, 1)

判断其边界条件可以使用 left_open 或者 right_open 来做：

>>> I.left_open
False
>>> I.right_open
False

对于两个集合之间的操作，可以参考以下方法：

I = sympy.Interval(0,1)
K = sympy.Interval(0.5,2)
>>> I.intersect(K)
Interval(0.500000000000000, 1)
>>> I.union(K)
Interval(0, 2)
>>> I-K
Interval.Ropen(0, 0.500000000000000)
>>> K-I
Interval.Lopen(1, 2)
>>> I.symmetric_difference(K)
Union(Interval.Ropen(0, 0.500000000000000), Interval.Lopen(1, 2))

实数集 $\mathbb{R}$ 在 SymPy 中用 sympy.S.Reals 来表示，自然数使用 sympy.S.Naturals，非负整数用 sympy.S.Naturals0，整数用 sympy.S.Integers 来表示。补集的计算可以用减号，也可以使用 complement 函数。

>>> sympy.S.Reals
Reals
>>> sympy.S.Reals-I
Union(Interval.open(-oo, 0), Interval.open(1, oo))
>>> I.complement(sympy.S.Reals)
Union(Interval.open(-oo, 0), Interval.open(1, oo))
>>> sympy.S.Reals.complement(I)
EmptySet()
>>> I.complement(K)
Interval.Lopen(1, 2)
>>> I.complement(sympy.S.Reals)
Union(Interval.open(-oo, 0), Interval.open(1, oo))

SymPy 的逻辑工具

在逻辑运算中，我们可以使用 $A, B, C$ 来代表元素。&, |, ~, >> 分别表示 AND，OR，NOT，imply。而逻辑运算同样可以使用 sympy.simplify_logic 简化。

A, B, C = sympy.symbols("A B C")
>>> sympy.simplify_logic(A | (A & B))
A
>>> sympy.simplify_logic((A>>B) & (B>>A))
(A & B) | (~A & ~B)
>>> A>>B
Implies(A, B)

SymPy 的级数工具

SymPy 的级数工具有一部分放在具体数学（Concrete Mathematics）章节了。有的时候，我们希望计算某个级数是发散的，还是收敛的，就可以使用 is_convergence 函数。考虑最常见的级数：

$\sum_{n=1}^{\infty}\frac{1}{n} = +\infty$ ，

$\sum_{n=1}^{\infty}\frac{1}{n^{2}} = \frac{\pi^{2}}{6}$ 。

>>> n = sympy.Symbol("n", integer=True)
>>> sympy.Sum(1/n, (n,1,sympy.oo)).is_convergent()
False
>>> sympy.Sum(1/n**2, (n,1,sympy.oo)).is_convergent()
True

如果想计算出收敛级数的值，加上 doit() 函数即可；如果想计算有效数字，加上 evalf() 函数即可。

>>> sympy.Sum(1/n**2, (n,1,sympy.oo)).evalf()
1.64493406684823
>>> sympy.Sum(1/n**2, (n,1,sympy.oo)).doit()
pi**2/6
>>> sympy.Sum(1/n**3, (n,1,sympy.oo)).evalf()
1.20205690315959
>>> sympy.Sum(1/n**3, (n,1,sympy.oo)).doit()
zeta(3)

除了加法之外，SymPy 也支持连乘，其符号是 sympy.Product，考虑

$\prod_{n=1}^{+\infty}\frac{n}{n+1}$ ，

$\prod_{n=1}^{+\infty}\cos\left(\frac{\pi}{n}\right)$ 。

>>> sympy.Product(n/(n+1), (n,1,sympy.oo)).is_convergent()
False
>>> sympy.Product(sympy.cos(sympy.pi/n), (n, 1, sympy.oo)).is_convergent()
True

SymPy 的 ODE 工具

在常微分方程（Ordinary Differential Equation）中，最常见的就是解方程，而解方程主要是靠 dsolve 函数。例如想求解以下的常微分方程：

$df/dx + f(x) = 0$ ,

$d^{2}f/dx^{2} + f(x) = 0$ ，

$d^{3}f/dx^{3} + f(x) = 0$ ，

可以使用 dsolve 函数：

>>> f = sympy.Function('f')
>>> sympy.dsolve(sympy.Derivative(f(x),x) + f(x), f(x))
Eq(f(x), C1*exp(-x))
>>> sympy.dsolve(sympy.Derivative(f(x),x,2) + f(x), f(x))
Eq(f(x), C1*sin(x) + C2*cos(x))
>>> sympy.dsolve(sympy.Derivative(f(x),x,3) + f(x), f(x))
Eq(f(x), C3*exp(-x) + (C1*sin(sqrt(3)*x/2) + C2*cos(sqrt(3)*x/2))*sqrt(exp(x)))

而常微分方程对于不同的方程类型也有着不同的解法，可以使用 classify_ode 来判断常微分方程的类型：

>>> sympy.classify_ode(sympy.Derivative(f(x),x) + f(x), f(x))
('separable', '1st_exact', '1st_linear', 'almost_linear', '1st_power_series', 'lie_group', 'nth_linear_constant_coeff_homogeneous', 'separable_Integral', '1st_exact_Integral', '1st_linear_Integral', 'almost_linear_Integral')
>>> sympy.classify_ode(sympy.Derivative(f(x),x,2) + f(x), f(x))
('nth_linear_constant_coeff_homogeneous', '2nd_power_series_ordinary')
>>> sympy.classify_ode(sympy.Derivative(f(x),x,3) + f(x), f(x))
('nth_linear_constant_coeff_homogeneous',)

SymPy 的 PDE 工具

在偏微分方程（Partitial Differential Equation）中，同样可以直接求解和判断偏微分方程的类型，分别使用函数 pdsolve() 和 classify_pde()。假设 $f = f(x,y)$ 是一个二元函数，分别满足以下偏微分方程：

$\partial f/\partial x + \partial f/\partial y =0$ ，

$\partial f/\partial x + \partial f/\partial y + f = 0$ ，

$\partial f/\partial x + \partial f/\partial y + f + 10 = 0$ 。

>>> f = sympy.Function("f")(x,y)
>>> sympy.pdsolve(sympy.Derivative(f,x)+sympy.Derivative(f,y),f)
Eq(f(x, y), F(x - y))
>>> sympy.pdsolve(f.diff(x)+f.diff(y)+f,f)
Eq(f(x, y), F(x - y)*exp(-x/2 - y/2))
>>> sympy.pdsolve(f.diff(x)+f.diff(y)+f+10,f)
Eq(f(x, y), F(x - y)*exp(-x/2 - y/2) - 10)

查看类型就用 classify_pde() 函数：

>>> sympy.classify_pde(f.diff(x)+f.diff(y)+f)
('1st_linear_constant_coeff_homogeneous',)
>>> sympy.classify_pde(f.diff(x)+f.diff(y)+f+10,f)
('1st_linear_constant_coeff', '1st_linear_constant_coeff_Integral')
>>> sympy.classify_pde(f.diff(x)+f.diff(y)+f+10,f)
('1st_linear_constant_coeff', '1st_linear_constant_coeff_Integral')

不过目前的 PDE 解法貌似只支持一阶偏导数，二阶或者以上的偏导数就不支持了。

SymPy 的数论工具

在数论中，素数就是一个最基本的概念之一。而素数的批量计算，比较快的方法就是筛法（sieve method）。在 sympy 中，同样有 sympy.sieve 这个工具，用于计算素数。如果想输出前100个素数，那么

>>> sympy.sieve._reset()
>>> sympy.sieve.extend_to_no(100)
>>> sympy.sieve._list
array('l', [2, 3, 5, 7, 11, 13, 17, 19, 23, 29, 31, 37, 41, 43, 47, 53, 59, 61, 67, 71, 73, 79, 83, 89, 97, 101, 103, 107, 109, 113, 127, 131, 137, 139, 149, 151, 157, 163, 167, 173, 179, 181, 191, 193, 197, 199, 211, 223, 227, 229, 233, 239, 241, 251, 257, 263, 269, 271, 277, 281, 283, 293, 307, 311, 313, 317, 331, 337, 347, 349, 353, 359, 367, 373, 379, 383, 389, 397, 401, 409, 419, 421, 431, 433, 439, 443, 449, 457, 461, 463, 467, 479, 487, 491, 499, 503, 509, 521, 523, 541, 547, 557, 563, 569, 571, 577, 587, 593, 599, 601, 607, 613, 617, 619, 631])

如果想输出一个区间内的所有素数，可以使用 primerange(a,b) 函数：

>>> [i for i in sympy.sieve.primerange(10,100)]
[11, 13, 17, 19, 23, 29, 31, 37, 41, 43, 47, 53, 59, 61, 67, 71, 73, 79, 83, 89, 97]

search() 函数是为了计算某个数附近是第几个素数：

>>> sympy.sieve.search(10)
(4, 5)
>>> sympy.sieve.search(11)
(5, 5)

如果只想获得第 n 个素数，则使用函数 sympy.ntheory.generate.prime(n) 即可。如果是希望计算 x 后面的下一个素数，使用 sympy.ntheory.generate.nextprime(x) 即可。判断 x 是否是素数，可以使用 sympy.ntheory.generate.isprime(x)。

>>> sympy.ntheory.generate.prime(10)
29
>>> sympy.ntheory.generate.nextprime(10)
11
>>> sympy.ntheory.generate.nextprime(11)
13
>>> sympy.ntheory.generate.isprime(11)
True
>>> sympy.ntheory.generate.isprime(12)
False

除此之外，SymPy 的数论方法还有很多，需要读者根据 SymPy 的官方文档自行探索。

SymPy 的范畴论工具

SymPy 还支持范畴论（Category Theory）的一些计算方法，在这里简要地列举一下。

>>> A = sympy.categories.Object("A")
>>> B = sympy.categories.Object("B")
>>> f = sympy.categories.NamedMorphism(A,B,"f")
>>> f.domain
Object("A")
>>> f.codomain
Object("B")

由于范畴论是数学的“黑话”，因此其余方法留给范畴论的科研人员自行翻阅。

总结：

整体来看，SymPy 是一个非常卓越的 Python 开源符号计算库。在符号计算领域，不仅支持常见的微积分，线性代数，几何运算，还支持集合论，微分方程，数论等诸多数学方向。后续笔者将会持续跟进并研究这一卓越的开源工具库。

参考文献：

Meurer A, Smith C P, Paprocki M, et al. SymPy: symbolic computing in Python[J]. PeerJ Computer Science, 2017, 3: e103.
GitHub：https://github.com/sympy/sympy
SymPy：https://www.sympy.org/en/index.html
Sympy 维基百科：https://en.wikipedia.org/wiki/SymPy
GreatX’s Blog：数值 Python：符号计算：https://vlight.me/2018/04/01/Numerical-Python-Symbolic-Computing/
SymPy 符号计算-让Python帮我们推公式：https://zhuanlan.zhihu.com/p/83822118
Python 科学计算利器—SymPy库：https://www.jianshu.com/p/339c91ae9f41

National University of Singapore

新加坡国立大学的数据科学与机器学习项目介绍

November 22, 2019 zr9558 1 Comment

新加坡国立大学（National University of Singapore）是一所综合性的大学，根据泰晤士报和世界大学排名来看，NUS 在整个亚洲的排名是非常靠前的。同时，数学系（Department of Mathematics）则是在理学院（Science）下的一个院系。

NUS数学系的介绍

新加坡国立大学数学系的前身可以追溯到 1929 年的 Raffles College。当时理学院开设了数学，化学，物理三门课程，不过总共也就只有十个学生和三位教师，其中有一位是数学教师。第一届数学系的领导（从 1931 年到 1959 年）是Alexander Oppenheim教授，他是在美国芝加哥大学获得的博士学位。从 1929 年开始，在新加坡的教育系统中，数学教育事业得到了巨大的发展，对现在的新加坡国立大学和南洋理工大学的建立起到了至关重要的作用。

NUS数学系首页

新加坡国立大学的数学系与国内的数学系有所不同。一般情况下，国内的数学系能够提供的专业包括数学与应用数学（Mathematics and Applied Mathematics），信息与计算科学（Information and Computing Science）与统计学（Statistics），有的时候会加上金融数学（Financial Mathematics）这一方向。而新加坡国立大学的数学系（Department of Mathematics）与统计系（Statistics）是分开的两个院系，虽然学生可以互相之间选择对方的课程，但是两者却是分属不同的院系。

NUS数学系本科专业

从数学系的首页来看，对于本科生而言，通常都有机会进行双学位的选择，例如：

计算机科学与数学；
经济学与数学；

因此，有不少的本科生都会有两个专业的学位证书。目标是为了让学生能够在未来从事数学研究和工业界的工作都打下坚实的基础。

NUS数学系硕士专业

NUS 的项目分成两块，coursework 项目和 Research 项目。第一个主要是以授课型的研究生为主，后者主要是为博士生或者准备攻读博士的人而准备的。对于硕士生而言，一般一两年就可以硕士毕业，只要修课的学分满了即可。

对于众多硕士生而言，进入 NUS 之前就要根据自身情况来选择一个合适的方向进行申请。对于硕士生的项目，数学系可以提供数学，金融数学等方向的课程，并且近期也提供了数据科学与机器学习专业（Data Science and Machine Learning）的课程。

项目介绍

申请条件

数据科学与机器学习项目是数学系，统计系，计算机系联合举办的为期一至两年的硕士生项目。期望学生的本科背景是数学方向，应用数学方向，统计与物理方向等。同时对学生的英语能力有一定的要求，希望对于母语不是英语的学生能够达到托福 85 分以上或者雅思 6.0 分以上的成绩。

项目课程要求

在课程安排方面，这个项目会有 20 个学分的课程，5 门核心课程，包括：

Introduction to Big Data for Industry;
Optimisation for Large-Scale Data-Driven Inference;
Foundations of Machine Learning/Theory and Algorithms for Machine Learning;
Cloud Computing;
DSML Industry Consulting and Applications Project.

而选修课方面包括机器学习，数据挖掘，大数据，计算机视觉，金融数学等方向的专业课程。同样也会要求选择 20 个学分的选修课程即可。

在学费方面，该项目是没有奖学金的，而一年的学费是 45000 新币，大约为 23 万人民币。

NUS申请流程

如果对这个项目感兴趣的同学，需要在 2020 年的 3 月 15 日之前提交申请，才有机会在 2020 年 8 月入学。

NUS数学系联系方式

如果有任何问题的话，可以考虑发邮件给 askmathpg@nus.edu.sg 或者拨打上面的电话号码。

PHD的生涯

读博的意义

November 4, 2019 zr9558 Leave a comment

读博的意义

在本科数学系，身边的学生大致分为几种情况：第一种学生是高中数学学的比较好，也没有什么特别爱好的专业，于是通过高考或者保送进入了数学系；第二种学生是为了转行计算机或者金融，就来数学系过渡一下，未来从事其他行业的工作；第三种学生是想读其他的专业，但是由于高考分数不足，被分配到了数学系；第四种学生就是为了攀登科学的高峰，立志要在数学界从事科研工作，在填报志愿的时候就果断进入数学系。一般情况下，学生在拿到本科学位之后，如果还想继续从事数学科研的话，就会选择继续攻读数学方向的博士学位，无论在国内还是国外继续学业。 640-6

在数学系本科或者硕士的时候，学生们通常都是按部就班的上课，做作业，上讨论班，或者每周集体分享一些论文或者学术资料。到了最后一个学期需要写论文的时候其实也不需要有太多的创新工作（更直接一点的话就是做不出来什么创新的工作），学生们大概就是整理一篇前沿的学术论文或者学术书籍就算完成了论文。在这两个阶段，导师基本上都会告诉学生，什么课题可以做得出来，什么课题比较难被做出来。或者就是导师其实已经知道某个问题的解决步骤和方法了，只需要学生执行完成就可以了。虽然在数学系的本科或者硕士都是按部就班的学习，但是这种学习是非常有必要的，毕竟不积跬步无以至千里，不积小流无以成江海。

QQ20191027-154338@2x

但是到了博士阶段就不一样了，导师给博士生的题目基本上就是学术界未知的问题。如果这个问题被解决了，那就不能够作为一个博士论文的题目。在某些极端情况下，导师连这个题目最终能不能够被学生做出来都不知道，只能够大致判断出这个问题的难易程度。数学系的纯数学方向和其他专业还不一样，某些专业跑跑数据改改模型就能够发表所谓的优秀论文了。而想要做纯数学的研究不仅需要通过天马行空的想象能力来寻找灵感，还需要有强大的逻辑推理能力把之前的想象整理成严谨的学术论文。而这种能力，仅仅靠本科或者硕士的训练是远远不够的，按部就班的读书考试也许能够应付小学，中学，大学，甚至研究生阶段，但是却无法完成博士阶段的任务目标。

QQ20191027-154259@2x

当年笔者在博士期间的研究方向是“复平面上多项式的 Julia 集合正测度”（PS：估计不是相关方向的人都不太读得懂这个题目的意思）。首先要知道复平面，其次要知道 Julia 集合，再次要知道正测度是什么意思。即使这些名词和概念都知道了，也不足以撬动这个题目。要想解决这个题目，除了导师的必要指导之外，还要自行去查阅各种论文资料来详细阅读。记得当年读过好几篇三五十页的 Annals of Mathematics（数学年刊）的论文，还读过一篇上百页的预印稿（关键是这篇论文有一个核心步骤还是错误的）。在此情况下，导师也只能够给博士生圈定一个论文的范围，告诉学生可以去参考其中的思路和方法，至于学生能不能够读懂这些文献，是导师无法保证的。而纯数学论文与其他专业最明显的不同就是满篇都是数学公式的推导，有很多地方充斥着“显然”，“易得”，“显而易见”诸多词汇。有一种情况是写论文的作者没有想明白，然后糊弄了一把；另外一种情况是这个地方真的是显而易见的，只是读者没有明白。纯数学论文的其中一页读上一两天并不是一件罕见的事情，一篇论文读一个学期能读明白也算完成了一件还不错的任务。PS：数学系有的教材就可以让学生一页读上好几天，例如 GTM 52。某些专业做不出来实验还可以把原因归结为材料不足，经费不足，设备不够，但是纯数学专业看不懂论文只能够把绝大部分原因归结为自己，因为草稿纸，笔，打印机，网络都是买得起的，客观上并不存在任何阻碍。

Julia_0.4_0.6 Julia-set_z2+c_-0.742_0.1_5000

在整个科研的过程中，一般都会经历：高峰->低谷->逐步爬坡回升->重新达到高峰的这样一个流程状态。正如：邓宁-克鲁格心理效应。从博士生的五年生涯来看，通常都会经历这几个不同的阶段：

1. 愚昧山峰：刚入学第一年的时候；

2. 低谷：第二年，第三年的时候；

3. 回升：第四年，第五年的时候。

而这种状态的经历是在读本科，硕士阶段都无法获得的，只有在经历了博士的求学之路才能够逐渐体会。当然，也有的导师也会帮学生写论文，整理论文，甚至送论文给学生，这样的学生可以在学生圈谋求一份教职，但是也永远体会不到这些心路历程了。这份心路历程在人生中是一场宝贵的财富，教会了如何在绝望之谷的时候迅速反弹回升，如何在逆境中持续成长，如何在被打击之后迅速恢复。

640-7

后来由于种种原因笔者没有继续从事博士期间的科研项目，到了企业从事普通的研发工作。在工作中，通常领导们也不会直接告诉下属项目的细节和执行方案，只会大致说出项目方向，该达成的目标是什么，下属们需要自行思考这个项目的背景是什么，为什么要启动这个项目，该怎么实现这个项目，怎么样做才能够把项目做好做完。这个过程跟读博期间的过程有点类似，但是又有所不同。相同的是都只能够拿到一个大致的方向，同时需要发挥自主能动性去做完细节，争取身边的各种资源，并且保证论文或者项目的落地。不同的是，工作期间并不会给下属很长的空窗期，并不会像博士期间给学生大几个月甚至一年的空闲时间去研究某个课题或者方向。基本上都是一旦制定了方向，就要立刻给出大致的方案和解决思路，然后给出项目的排期和执行计划，并且给出项目的风险点是什么，最终开始逐步落实方案和具体的工作内容。不过，在工作中，在大多数情况下也不会要求一定要做一个行业第一的产品，而且很多项目由于种种原因也无法做到行业第一。有的时候，行业第一的量化标准还不太好做，如果有明确规则的，例如挣钱的多少，日活月活多少，比赛的第一，这种是可以制定出来的。但是大多数情况下其实也不太好制定出来。不过在工作中还是可以精益求精，以行业的第一为目标来逐步提升的。可以在现有的情况和条件下把项目做得尽量完善，在提升效率方面比之前能够更加好，在节省人力方面能够做到较少人力就维护较多的事情。

在工作之后，虽然也很难体会到一篇论文读半年的时候，但是在博士期间学会的资料整理，收集，汇总，甚至日常的笔记和资料总结，都是在工作中非常常用的技能。甚至一些时间管理和项目管理的经验，也是在工作中必备的技能之一。

640-2

640-4

很难说读博这件事情对所有的人都有意义，因为每个人都有着不同的生长环境和家庭条件，在如此复杂的情况下，给出一个固定的答案其实没有任何意义。何况每个人的想法也会随着时间的迁移，环境的变化而做出改变。之前想做的一些事情可能后面就不感兴趣了，一直想追求的一些目标可能也不再继续追求了。只能够根据具体的情况，在不同的阶段来定制化地给每个人相应的建议和支持（个性化推荐）。无论做什么样的事情，最好都要保持一种乐观向上的心态，不能够一直处于一种低估和不能自拔的阶段。如果在本科，硕士，甚至博士期间都没有找好自己的定位，也可以尝试在工作之后逐步把自己拉上正规，在期间调整自己的状态和节奏，让自己的幸福感逐步提升。无论在哪里，无论最终是否做科研，最终都是为了过上更好的生活。

张戎

2019年10月27日

时间序列

深度学习在时间序列分类中的应用

September 6, 2019 zr9558 Leave a comment

本篇博客将会分享几篇文章，其内容主要集中在深度学习算法在时间序列分类中的应用。

无论是图像分类，文本分类，还是推荐系统的物品分类，都是机器学习中的常见问题和应用场景。同样的，时间序列的分类问题也是研究时间序列领域的重要问题之一。近期，神经网络算法被用于物体识别，人脸识别，语音分类等方向中，于是有学者用深度学习来做时间序列的分类。

假设

$X=\{x_{1},\cdots,x_{n}\}$

是一个长度为 $n$ 的时间序列，高维时间序列

$X = \{X^{1},\cdots,X^{M}\}$

则是由 $M$ 个不同的单维时间序列而组成的，对于每一个 $1\leq i\leq M$ 而言，时间序列 $X^{i}$ 的长度都是 $n.$ 而时间序列的分类数据通常来说都是这种格式：数据集

$D=\{(X_{1},Y_{1}),(X_{2},Y_{2}),\cdots,(X_{N},Y_{N})\}$

表示时间序列与之相应的标签，而 $Y_{i}$ 是 one hot 编码，长度为 $K$ （表示有 $K$ 个类别）。

整体来看，时间序列分类的深度学习方案大体是这个样子的：输入的是时间序列，通过某个神经网络算法进行端到端的训练，最后输出相应的分类概率。

而做时间序列分类的深度学习算法分成生成式（Generative）和判别式（Discriminative）两种方法。在生成式里面包括 Auto Encoder 和 Echo State Networks 等算法，在判别式里面，包括时间序列的特征工程和各种有监督算法，还有端到端的深度学习方法。在端到端的深度学习方法里面，包括前馈神经网络，卷积神经网络，或者其余混合模型等常见算法。

深度学习算法在时间序列分类中的应用：Baseline

这一部分将会介绍用神经网络算法来做时间序列分类的 Baseline，其中包括三种算法，分别是多层感知机（MLP），FCN（Fully Convolutional Network）和 ResNet。其论文的全名是《Time Series Classification from Scratch with Deep Neural Networks: A Strong Baseline》。这篇论文中使用的神经网络框架如下图所示：

DNN_Baseline_结构1.png

多层感知机（MLP）模型使用了全连接层，每个隐藏层大约 500 个神经元，然后使用 ReLU 作为激活函数，同时使用 Dropout 来防止过拟合，最后一层是 Softmax 层。MLP 中一个基础的块包括：

$\tilde{x} = f_{dropout, p}(x)$ ,

$y = W\cdot \tilde{x} + b$ ,

$h = ReLU(y)$ .

除了前馈神经网络之外，全卷积网络（FCN）同样可以作为时间序列的特征提取工具，一个卷积块包括：

$y = W \otimes x + b$ ,

$s = BN(y)$ ,

$h = ReLU(s)$ ,

在这里， $\otimes$ 指的是卷积算子，BN 指的是 Batch Normalization，ReLU 则是激活函数。

Residual Network 是在 FCN 的基础上进行的改造。令 $Block_{k}$ 来表示第 $k$ 个卷积块，而 Residual 块就定义为：

$h_{1} = Block_{k_{1}}(x)$ ,

$h_{2} = Block_{k_{2}}(h_{1})$ ,

$h_{3} = Block_{k_{3}}(h_{2})$ ,

$y=h_{3} + x$ ,

$\hat{h}=ReLU(y)$ .

其中， $k_{1} = 64, k_{2} = 128, k_{3} = 128$ 。

评价指标

Mean Per Class Error (in Multi-class Classification only) is the average of the errors of each class in your multi-class data set. This metric speaks toward misclassification of the data across the classes. The lower this metric, the better.

模型的评价指标使用的是 Mean Per-Class Error，指的是在多分类场景下，每一类（Class）错误率的平均值。换句话说，一个数据集 $D=\{d_{k}\}_{1\leq k\leq K}$ 是由 $K$ 个类的元素构成的，每个类的标签是 $C=\{c_{k}\}_{1\leq k\leq K}$ ，通过模型其实可以计算出模型对每一个类的错误率 $e_{k}$ ，那么模型的 MPCE 就是： $MPCE= \sum_{1\leq k\leq K} e_{k}/K$ .

其实验结论是：

DNN_Baseline_实验数据1.png

MSCNN

MSCNN 的全称是 Multi-Scale Convolutional Neural Networks，相应的论文是《Multi-Scale Convolutional Neural Networks for Time Series Classification》。

在时间序列的分类算法里面，通常来说，可以分成以下几种：

基于距离的方法（distance-based methods）：kNN，SVM（相似核），DTW；
基于特征的方法（feature-based methods）：SVM，逻辑回归等；
基于神经网络的方法（neural network-based methods）：CNN 等；

正如前文所提到的，一条时间序列通常可以写作 $T=\{t_{1},\cdots,t_{n}\}$ ，其中 $t_{i}$ 表示在时间戳 $i$ 下的取值，并且时间序列 $T$ 的长度是 $n$ 。在时间序列分类的场景下，每一条时间序列对应着唯一的一个标签（label），也就是说 $D=\{(T_{i},y_{i})\}_{i=1}^{N}$ 。其中 $D$ 集合里面包含 $N$ 条时间序列，每条时间序列 $T_{i}$ 对应着一个标签 $y_{i}$ 。 $y_{i}$ 表示分类值集合 $\mathcal{C} = \{1,\cdots,C\}$ 中的元素， $C\in \mathbb{Z}^{+}$ 。

MSCNN 的整体结构：

在 Multi-Scale Convolutional Neural Network（MSCNN）中，包括几个串行的阶段，

变换阶段（Transformation Stage）：包括恒等变换，下采样，谱变换等变换方式，每一种方式都是一个分支，并且也是卷积神经网络的输入；
局部卷积（Local Convolution Stage）：使用卷积层来对不同的输入提取特征，不同的输入分支之间是相互独立的，输出的时候都会经过一个最大值池化（max pooling）的过程；
整体卷积（Full Convolution Stage）：把上一步提取到的特征进行拼接（concatenate），然后使用全连接层并且加上一个 softmax 层来做多分类。

如下图所示，MSCNN 是一个端到端的训练网络结构，所有参数都是通过后向传播算法得到的。

MSCNN结构1.png

首先来看神经网络的第一步，变换阶段（Transformation Stage），也就是神经网络的多尺度的输入。在不同的尺度下，神经网络能够提炼到不同类型的特征。长期的特征（long-term features）反映了时间序列的整体趋势，短期的特征（short-term features）反映了时间序列的局部的微妙变化。要想判断时间序列的形状，不仅要参考整体的特征，也要参考局部的特征，这两者对于判断时间序列的形状都具有一定的辅助作用。

在 Transformation Stage，identity map 指的是恒等变换，也就是说时间序列是原封不动的作为神经网络的输入数据。对于 Smoothing Transformation，指的就是对时间序列进行必要的平滑操作，将新的时间序列作为神经网络的输入数据。在这种情况下，我们可以对时间序列 $T=\{t_{1},\cdots,t_{n}\}$ 进行移动平滑，i.e.

$T^{\ell}=(x_{i}+x_{i+1}+\cdots+x_{i+\ell-1})/\ell, 0\leq i\leq n-\ell+1$ ,

其中的 $\ell\in \mathbb{Z}^{+}$ 表示窗口长度。对于不同的窗口长度 $\ell$ ，我们可以的到不同的时间序列平滑序列，但是它们的长度都是一样的，都是原始的时间序列长度 $n$ 。

而下采样（down sampling）指的则是对时间序列的间隔进行抽样操作。假设时间序列 $T=\{t_{1},\cdots,t_{n}\}$ ，下采样的比例是 $k$ ，也就是说我们每隔 $k$ 个点保留时间序列的取值，i.e.

$T^{k} = \{t_{1+k\cdot i}\}, 0\leq i \leq [(n-1)/k]$ .

用这种方法，我们可以对 $k=1,2,3,\cdots$ 来进行下采样的时间序列提取。在进行了恒等变换，平滑变换，下采样之后，时间序列就可以变成多种形式，作为神经网络的输入。

其次，在神经网络部分，本文使用了一维（1-D）的卷积层和最大值池化的方法来提取特征，并且在局部卷积阶段之后把提炼到的抽象特征进行拼接（concatenate）。拼接完了之后，持续使用卷积层和池化层进行特征的提取，然后使用全连接层（fully connected layers）和 softmax 层来进行时间序列类别的预测。

数据增强

在深度学习里面，由于是端到端的训练网络，因此是需要相对多的样本数据的，于是有的时候需要进行数据增强（data augmentation）。也就是在现有的基础上获得更多的训练数据。对于时间序列 $T=\{t_{1},\cdots,t_{n}\}$ ，可以定义一个子序列：

$S_{i:j} = \{t_{i}, t_{i+1},\cdots,t_{j}\}, 1\leq i,\leq j\leq n$ ，

对于正整数 $s\in \mathbb{Z}^{+}$ ，可以生成 $n-s+1$ 个子序列如下所示：

$Slicing(T,s) = \{S_{1:s},S_{2:s+1},\cdots,S_{n-s+1:n}\}$ ,

这些子序列的标签与原始的时间序列 $T$ 是一样的。

本文用到的数据集情况如下表所示：

MSCNN数据集1.png

实验数据如下图所示：

MSCNN实验1 MSCNN数据集2

结论

本文使用了 MCNN 来对变换之后的时间序列进行特征提取，并且进行了端到端的模型训练。并且也讨论了卷积神经网络使用在 shapelet learning 上的一些逻辑和方法，然后解释了 MSCNN 在时间序列分类上能够有不错表现的原因。但是所有的 TSC 数据集都不算特别大，对端到端的训练模式有一定的限制。

GASF 和 GADF 方法

这篇文章《Imaging Time Series to Improve Classification and Imputation》介绍了如何把时间序列转换成图像，包括 GASF 方法和 GADF 方法。

假设时间序列是 $X = \{x_{1},\cdots, x_{n}\}$ ，长度是 $n$ ，我们可以使用归一化方法把时间序列压缩到 $[0,1]$ 或者 $[-1,1]$ ：

$\tilde{x}_{0}^{i} = (x_{i} - \min(X))/(\max(X) -\min(X))$ ,

$\tilde{x}_{-1}^{i} = ((x_{i} - \max(X)) + (x_{i}-\min(X)))/(\max(X) - \min(X))$ ,

此时的 $\tilde{x}_{0}^{i}\in[0,1], \forall 1\leq i\leq n$ ， $\tilde{x}_{-1}^{i} \in [-1,1],\forall 1\leq i\leq n$ 。于是可以使用三角函数来代替归一化之后的值。下面通用 $\tilde{x}_{i}$ 来表示归一化之后的时间序列，令 $\phi_{i} = \arccos(\tilde{x}_{i})$ ， $\tilde{x}_{i} \in [-1,1]$ ， $1\leq i\leq n$ 。因此， $\phi_{i}\in[0,\pi]$ ，于是， $\sin(\phi_{i}) \geq 0$ 。

定义矩阵 GASF（Gramian Angular Summation Field） 为

$GASF = [\cos(\phi_{i}+\phi_{j})]_{1\leq i,j\leq n}$

于是，

$GASF = [\cos(\phi_{i})\cdot \cos(\phi_{j}) - \sin(\phi_{i})\cdot \sin(\phi_{j})]_{n\times n}$

令 $\tilde{X}=(\cos(\phi_{1}),\cdots,\cos(\phi_{n}))^{T}$ ，可以得到

$GASF = \tilde{X} \cdot \tilde{X}^{T} - \sqrt{I - \tilde{X}^{2}} \cdot \sqrt{I - \tilde{X^{T}}^{2}}$

以上的都是 element 乘法和加法， $I$ 表示单位矩阵。它的对角矩阵是

$diag(GASF) = \{\cos(2\phi_{1}),\cdots, \cos(2\phi_{n})\}$

$= \{2\cos^{2}(\phi_{1})-1,\cdots,2\cos^{2}(\phi_{n})-1)\} = \{GASF_{ii}\}_{1\leq i\leq n}.$

如果是使用 min-max normalization 的话，是可以从 diag(GASF) 反推出 $\tilde{x}_{i}$ 的。因为， $2 \tilde{x}_{i}^{2} - 1 = 2\cos^{2}(\phi_{i}) - 1 = GASF_{ii}$ ，可以得到 $y_{i} = \sqrt{(GASF_{ii}+1)/2}$ 。

定义 GADF（Gramian Angular Difference Field）如下：

$GADF = [\sin(\phi_{i}-\phi_{j})]_{1\leq i,j\leq n}$

$= [\sin(\phi_{i})\cdot cos(\phi_{j}) - \cos(\phi_{i})\cdot\sin(\phi_{j})]_{1\leq i,j\leq n}$

$= \sqrt{1-X^{2}}\cdot X^{T} - X \cdot \sqrt{1-(X^{T})^{2}}$ .

Markov Transition Field（MTF）

除了 GSAF 和 GSDF 之外，《Imaging Time Series to Improve Classification and Imputation》，《Encoding Time Series as Images for Visual Inspection and Classification Using Tiled Convolutional Neural Networks》，《Encoding Temporal Markov Dynamics in Graph for Time Series Visualization》也提到了把时间序列转换成矩阵 Image 的算法 MTF。在 pyts 开源工具库里面，也提到了 MTF 算法的源码。

假设时间序列是 $X = \{x_{1},\cdots,x_{n}\}$ ，我们把它们的值域分成 $Q$ 个桶，那么每一个 $x_{i}$ 都可以被映射到一个相应的 $q_{j}$ 上。于是我们可以建立一个 $Q\times Q$ 的矩阵 $W$ ， $w_{ij}$ 表示在桶 $j$ 中的元素被在桶 $i$ 中的元素跟随的概率，也就是说 $w_{ij} = P(x_{t}\in q_{i}|x_{t-1}\in q_{j})$ ，同时，它也满足 $\sum_{j=1}^{Q}w_{ij} =1$ 。于是，得到矩阵 $W = (w_{ij})_{1\leq i,j\leq Q}$ 。

除此之外，我们也能够计算一个迁移概率矩阵 $M$ 。其中 $m_{ij}$ 表示桶 $i$ 中的元素迁移至桶 $j$ 中的概率 $P(q_{i}\rightarrow q_{j})$ ，同样有 $\sum_{1\leq j\leq Q} m_{ij} =1$ 。因此，我们同样可以构造出一个 $Q\times Q$ 的矩阵将时间序列可视化。

时间序列的降维方法有两种：

分段聚合（PAA）：使用局部平均等方法，把时间序列进行降维；
核变换（Kernel）：使用 Bivariate Gaussian 核或者均值核来把时间序列进行降维。

在把时间序列进行可视化之后，对于时间序列分类的场景，就可以使用 CNN 的技术方案来做了。如下图所示：

其实验数据效果如下：

tildeCNN实验数据1.png

Time Le-Net

在本篇文章《Data Augmentation for Time Series Classification using Convolutional Neural Networks》中，主要用到了卷积神经网络来做时间序列的分类。

除此之外，也使用了不少数据增强（Data Augmentation）的技术。包括前面提到的 Window Slicing（WS）方法。也考虑了 Warping 的变换技巧，例如 Warping Ratio = 1/2 或者 2。这种时间扭曲指标比率可以通过交叉验证来选择。该方法叫做 Window Warping（WW）技术。

另外也有其余论文使用卷积神经网络做时间序列分类，例如《Convolutional neural networks for time series classification》，如下图所示：

Multi-Channels Deep Convolutional Neural Networks

在高维时间序列的分类中，有人提出用多通道的卷积神经网络来进行建模。

整体来看，分成四个部分。前三个部分作为特征提取工具，最后一层作为分类工具。

Filter Layer：
Activation Layer：
Pooling Layer：
Fully-Connected Layer：

实验对比数据如下：

结论：

在本篇博客中，列举了一些深度学习算法在时间序列分类中的应用，也介绍了部分数据增强的方法和时间序列数据变换的方法。从以上各篇文章的介绍来看，深度学习在时间序列分类领域上应该是大有可为的。

参考资料：

Wang Z, Yan W, Oates T. Time series classification from scratch with deep neural networks: A strong baseline[C]//2017 international joint conference on neural networks (IJCNN). IEEE, 2017: 1578-1585.
Cui Z, Chen W, Chen Y. Multi-scale convolutional neural networks for time series classification[J]. arXiv preprint arXiv:1603.06995, 2016.
Wang Z, Oates T. Imaging time-series to improve classification and imputation[C]//Twenty-Fourth International Joint Conference on Artificial Intelligence. 2015.
Wang Z, Oates T. Encoding time series as images for visual inspection and classification using tiled convolutional neural networks[C]//Workshops at the Twenty-Ninth AAAI Conference on Artificial Intelligence. 2015.
Liu L. Encoding Temporal Markov Dynamics in Graph for Time Series Visualization, Arxiv, 2016.
Fawaz H I, Forestier G, Weber J, et al. Deep learning for time series classification: a review[J]. Data Mining and Knowledge Discovery, 2019, 33(4): 917-963.
Zhao B, Lu H, Chen S, et al. Convolutional neural networks for time series classification[J]. Journal of Systems Engineering and Electronics, 2017, 28(1): 162-169.
Le Guennec A, Malinowski S, Tavenard R. Data augmentation for time series classification using convolutional neural networks[C]. 2016.

时间序列

时间序列的标签

August 16, 2019 zr9558 Leave a comment

本篇文章是为了介绍一种基于少量样本标记而获得更多样本的方法，论文的原文是《Label-Less: A Semi-Automatic Labeling Tool for KPI Anomalies》，是清华大学与多家公司（必示科技，中国建设银行等）的合作论文。

在时间序列异常检测中，因为标注的成本比较大，于是需要寻找一种较少而高效地标注时间序列异常点的方法。在该论文中，Alibaba，Tencent，Baidu，eBay，Sogou提供了上千条时间序列（每条时间序列大约是2-6个月的时间跨度），作者们进行了 30 条 KPIs 的标注工作。但是其标注成本依旧是很大的，于是作者们想到了一种异常相似搜索（anomaly similarity search）的算法，目标是对已经标注好的时间序列异常模式进行模版搜索。目的就是达到 label-less，也就是较少的标注而获得更多的标注数据。

在本篇论文中，在异常检测的过程中，作者们使用了时间序列的预测模型（time series prediction models）来获得时间序列的特征，使用了孤立森林（Isolation Forest） 来对时间序列的特征来做无监督的异常检测。并且其效果由于 one class svm 算法和 local outlier factor 算法。在搜索的部分，作者使用了加速版的 DTW 算法（accelerated dynamic time warping approach）来做相似度的搜索和模式的匹配。其中也尝试了各种技巧和方法，包括 constrained DTW，LB Keogh 方法，early stopping 算法等工具。

整个 Label-Less 的架构图如下表示：

其中的 Operators 指的是业务运维人员，面对着无标记的多条时间序列曲线。系统首先会进行无监督的异常检测算法啊，包括时间序列的预处理（归一化等）操作，然后使用差分（Difference），移动平均算法（moving average），带权重的移动平均算法（weighted moving average），指数移动平均（ewma），holt winters，ARIMA 等算法来做特征的提取。此时，对于不同的时间序列预测工具，我们可以得到不同的预测值，然后把预测值减去实际值并且取绝对值，就得到时间序列的误差序列。i.e. $|p_{i} - x_{i}|$ 就作为数据点 $x_{i}$ 的特征。

在这种情况下，由于用了六个时间序列预测算法，因此原始的时间序列 $X (n\times 1)$ 就可以变成特征矩阵 $X' (n\times 1)$ 。对于特征矩阵 $X'$ 可以使用 isolation forest 来做无监督的异常检测并且做阈值的设定；如下图所示：

而另外的一部分的异常相似搜索（anomaly similarity search）是在第一部分的基础上在做的，Unsupervised Anomaly Detection 会输出疑似异常或者候选异常，并且基于已知的异常模板（Anomaly Template）进行相似度的匹配，此时可以使用 accelerated DTW 算法，选择出最相似的 Top-K 异常，然后运维人员进行标注，得到更多的样本。

由于，对于两条长度分别是 $m$ 和 $n$ 的时间序列，DTW 相似度算法的时间复杂度是 $O(mn)$ ，因此在搜索的时候需要必要的加速工作。在这种地方，作者们使用了 LB-Kim，LB-Keogh，LB-Keogh-Reverse 算法来做搜索的加速工作。而这些的时间复杂度是 $O(m+n)$ 。整体的思路是，如果两条时间序列 $q$ 和 $c$ 的 LB-Kim，LB-Keogh，LB-Keogh-Reverse 的下界大于某个阈值，则不计算它们之间的 DTW 距离。否则就开始计算 DTW。并且在计算 DTW 的时候，如果大于下界，则会提前终止（early stopping），不会继续计算下去。如果都没有大于阈值，则把这个候选曲线和 dist 距离放入列表，最后根据列表中的 dist 来做距离的逆序排列。

整体流程如下：

AnomalySimilaritySearch

其运行速度也比直接使用 DTW 快不少：

Label-Less 的交互页面如下所示：

图（a）表示使用无监督算法获得的疑似异常；

图（b）表示使用异常搜索算法获得的异常结果。

下图则表示模板， $m$ 表示模板的长度， $c$ 表示相似的异常候选集个数；

Fig11

总结：

整体来看，本文提供了一种通过少量人工标注，无监督算法和相似度算法来获得更多样本的方法。在候选的时间序列条数足够多的时候，是可以进行时间序列的相似度匹配的。这给未来在运维领域提供海量的时间序列标注数据给予了一定的技术支持。

互联网的生涯

怀念 Jett

August 11, 2019 zr9558 Leave a comment

之前笔者写过不少文章，其中不仅包括技术和数学文章，还有一些日常学习工作的心得，甚至还写过几篇关于旅游的流水账，但是却很少写关于人物的文章。最近我打算写一篇关于人物的文章，用于纪念在工作中一起战斗过的伙伴 Jett。

许多年前，我还在学校读书的时候，就有前辈告诉我，交朋友的时候基本上就只能够在学校，离开了学校之后在职场中想交到合适的朋友是比较困难的。当时作为一枚职场菜鸟的我，也只能够表示对这句话的认同。直到在工作后遇到了 Jett，我对这个观点就有了自己的看法。

屏幕快照 2019-08-10 下午10.02.49

我是2015年7月作为校招生进入公司，Jett 则是2015年10月通过社招进公司的。当年虽然两人身处同一个部门，但是由于不在一个中心，也不是同一条业务线，并且这两个中心平时也没啥日常工作的交集，因此我和 Jett 同处一个部门许久却互相不认识对方。就好比在学校里面，虽然身处同一个院系，甚至同一栋宿舍楼，由于种种原因而不认识对方。

随着工作时间的增加，到 2016 年的时候笔者就已经过了新手保护期。所谓新手保护期，就是每一个新人其实都有一段适应公司的过程，公司也会给人一定的保护时间，但是一旦过了这个时间就不算新人了。大约到了2016年5月份的时候，部门内部启动了一个用机器学习的方法来做安全对抗的项目。当年的安全中心并没有机器学习的从业者，因此需要从其他中心抽调人手去支持他们的项目。在一些机缘巧合之下，我作为数据团队的代表，就和安全团队的同事们同在安全大数据项目组。在公司里面，机器学习的人员通常都是和开发人员通力合作，各司其职，共建一个项目。在项目初期，其实 Jett 也不负责与我对接，后续随着项目的进一步发展和迭代，Jett 就成为了安全中心与我对接的人。

虽然 Jett 是和我对接的人，但是 Jett 主要也是负责平台开发，之前 Jett 在进入公司之前好像没有做过前端，但是在项目初期好像做了一些前端方面的工作，在邮件里面得到了领导的赞扬。这个应该是我对 Jett 的第一印象了，那就是一个靠谱的开发人员。当时虽然是成立了项目组，但是还是作为一个试点的工作来做，无论是算法的调研和平台的工作，都由我们两个来完成。这些项目在初期都是一大堆人在里面，但是随着项目的进展，人数都会越来越少，然后领奖的时候又出现了一大堆人。在项目的初期，我当时应该也只会 SQL，连 Python 都不算很熟悉，更没有开发过大型项目了。即使是作为一个试点的项目，其实个人也有着巨大的压力，毕竟这个方向在公司当年还处于一个相对不那么清晰的状态，也不知道该怎么做，甚至也不知道有什么资料可以查。不过作为一个 PHD，虽然写代码的能力不行，但是读论文和搜索资料的能力还是有点的。当时记得花了很多时间在网上找来找去，总算找到一两家号称做这个方向的创业公司，运气够好的是当时它们居然都发表了相关的论文。后来把它们的论文整理成资料和PPT，跟相关的同事沟通了一下，就打算复现这几篇论文。

屏幕快照 2019-08-10 下午10.01.49

当然，做项目总是有着各种各样的风险，况且当年本人的技术实力确实也不怎么样，写过不少的 SQL，做过很多数据分析的工作，查过一些数据方面的问题，至于其他方面也就那么回事。在项目的调研阶段和初期，大约是5月份-8月份，其实也没有什么特别大的产出，主要也是在积攒技术经验，把 Python 的各个工具库熟悉了一下。不过在公司里面总是有考核的，到了9月份的时候压力就比较大了，毕竟最终还是要对上线的效果负责的。在9月底的时候，我和 Jett 被拉去星巴克喝了一杯咖啡，然后被下达了任务和命令，也就是在十月份完成机器学习的效果指标。

屏幕快照 2019-08-10 下午9.59.46

既然被拉到星巴克去喝咖啡，也就是说明在9月份的时候，模型的效果都不算太好，估计再这样干下去也没有办法完成既定的目标。不过目标这种事情制定的时候都是根据当时的情况来制定的，很有可能会随着时间的变化而产生变化。虽然说目标总是在发生变化，但是整体的大思路是没有改变的。当时，作为开发人员的 Jett 问我要不要修改 XGBoost 的模型源码，我连忙说不用不用，其实就算让我改，我也没有本事修改别人的源码。既然不会修改模型的源码，那就只能够从数据和特征来入手了。在机器学习领域，只要把数据和特征处理好了，基本上就能够保证模型的效果和质量。于是 Jett 和我当时把成千上万条数据一起看了一遍，结果是发现正负样本有一部分混在了一起，于是模型的效果无论如何都做不好。

屏幕快照 2019-08-10 下午10.05.45 — XGBoost 的 Github 主页

对于在校学生而言，通常来说都是 Python 跑完模型，然后得到一个模型文件，用它继续离线预测就可以得到最终的结果。但是在工业界，很多机器学习项目都需要进行上线的工作。这种时候只靠一个机器学习人员的战斗力是无法解决问题的，不可避免地需要有开发人力的介入。此时 Jett 发挥了作为一位开发人员的强大战斗力，一个人就能把线上的代码全部完成，没有让我撰写任何一行 C++ 的代码。在码农界有一种“结伴编程”的说法，也就是两个人共同搞一份代码，共同搞一个项目。我是负责数据处理和离线模块，Jett 是负责平台开发。有一次在核对数据的时候，把 Python 的预测结果和 C++ 的预测结果进行核对。也就是为了交流方便，那次我把 Mac Air 搬到了 Jett 的桌子上，两个人并排坐在一起，一起核对数据的准确性和可靠性。

不过核对数据只是机器学习项目的第一步，并且也是长期需要做的一步，因为数据总是会出现各种各样的问题。即使数据没有问题，也不代表最终的模型效果达标。到了2016年10月中旬的时候，我俩的压力也已经很大，每天晚上都在万利达五楼加班到十点以后。每天我都在训练模型和数据分析，Jett 每天也在做平台开发的工作，虽然我不懂他在做什么，但是总是感觉很忙的样子。而且项目经理每两天都会催一次进度，顺带着会有各种大大小小的会议。其实当时 Jett 是作为业务部门的人，压力更大的是压在他那一边。在思考了大半个月之后，在一个星期日的夜晚，我去华润万家买了一堆小本子，把整个项目的细节仔仔细细地思考了一遍，突然灵光一现，发现其实有一个优化点，于是就顺手写在了小本子上。到了周一的时候，我把小本子上面的想法去实现了一遍，结果效果瞬间提升，业务指标瞬间完成。当我把这个好消息告诉 Jett 的时候，Jett 还在 RTX 上说了一句：“等我去哭一会”。最终我俩猛加班两天，把模型往线上一扔，效果直线上升，而此刻的报表系统早已齐备，瞬间就看到了效果的提升，甩掉前面的模型五条大马路。

在2017年初的时候，模型的整体效果已经很好了，整体的路线其实已经基本走通。因此后续更多的是平台建设方面的工作，这一块做得更多一点的是 Jett，最终他还真是完成了无算法人员参与就能够自动接入各种业务的全自动流程。不过随着时间的流逝，我俩都各自都被抽调去做其它的项目，合作的时间和机会在项目成功之后就变得越来越少了。万幸的是，在2017年9月5日-2017年9月9日，我俩共同去了台湾进行了一次团建。参见《五日台湾行》。

640

后来，由于种种原因，我也不再继续从事安全方面的工作，与 Jett 的交集也越来越少。虽然后面也有吃饭聚餐等时候，但是却没有了在一个项目中一起为了一个目标奋斗的机会，心中不免有少许遗憾。随着去年的大调整，Jett 和我已经不在同一个部门。今年 Jett 由于家庭原因离开了公司，恐怕短期内已经再也没有通力协作，互帮互助的机会了。也许这个项目是我俩共同做的最后一个项目，但是当年在万利达五楼共同奋斗的时光，却是一段珍贵的回忆。

Nanjing University

高考志愿—谈一谈数学专业

July 1, 2019 zr9558 Leave a comment

最近，经常在知乎或者其他平台上看到

“高考填写志愿要不要填写数学专业？”
“读计算机专业是不是需要先去数学系？”
“从数学系转行到金融行业的前景怎么样？”

这一类的问题。最近也正值高考填志愿的时期，于是在这里说一下自己的一些想法。

其实这个问题很难给出一个标准的答案，因为每个人的个人条件，家庭条件，就业规划都完全不一样。推荐的结果应该正如推荐系统一样，根据每个人的情况来定制化，而不是给出一个结果供所有人参考。

一般情况下，每个人都有自己的发展方向或者职业规划，基本上进入数学系有这几个常见的原因：

想做数学科研，成为数学工作者；
为了转行做准备，准备研究生阶段进入金融或者计算机；
实在不知道想学啥，于是选择数学试试看。
其他原因。

Case 1:

第一种情况，想成为数学工作者的应届生其实只能选择数学系。目前的现代数学高度专业化，已经不是一个外行人，跨专业的人或者普通职业的人能够研究的了。本科所传授的数学知识大约也就是在20世纪中期而已，后面几十年的科研工作就很难写进本科教材了。因此，对于有志向成为数学家的学生而言，进入数学系是唯一的选择。当然其实也有大神从别的专业转过来从事数学工作的，并且也能够做得很好，但是对于绝大多数人而言是不现实的。

对于这一批少数的学生而言，最好的方式就是按照数学系的标准培养模式来做。也就是听课，刷题，复习，钻研，考试这条路，每天花十个小时学习，高强度的学习几年的本科数学。然后想办法进入顶尖高校的数学系，继续从事数学科研几年时间。最后进入顶尖高校的数学系或者科研站，从事数学研究。

Case 2:

第二种情况，经常都会遇到一个问题，那就是“是否有必要先读数学，然后在研究生阶段转行计算机或者金融？”先说结论：个人感觉是完全没有必要的。理由有以下几点：

金融或者计算机同样会开设各种各样的数学课，并且是结合专业的情况来开设相关课程；
数学系并不会专门提供金融或者计算机的课程；
数学系的思维方式和培养模式与金融或者计算机不太一样。

在这里，如果是为了打基础而先进入数学系再转行去其他专业其实是没有必要的。因为其他专业同样可以提供数学课，完全可以在其他专业学好数学课，没有必要专门跑来数学系学习四年数学。当然，由于高考填写志愿，或者为了进入一个更好的大学，只能够选择数学系而放弃其他专业其实是可以理解的。但是，在本科期间一定要补充其他专业技能，否则最后还是会比较吃亏。

除此之外，数学系很少会去开设金融或者计算机方面的课程，大部分开设的还是数学课，当然也不排除有的数学系会开设一大堆计算机系课程的情况。那么就来看一下数学系与计算机系的课程安排：

数学系的课程：

数学分析，高等代数，解析几何，C++，离散数学，常微分方程，偏微分方程，抽象代数，复变函数，实变函数，泛函分析，数值计算，偏微分方程数值解，拓扑学，微分几何，概率论与数理统计，随机过程等。

计算机系的课程：

微积分，线性代数，离散数学，数据结构与算法，数字电路，计算机组成原理，操作系统，编译原理，计算机网络，数据库原理，软件工程，汇编语言等。

从这两个课表的对比情况来看，如果要想从数学系转行到计算机系，那么基本上要把计算机的一些基础知识课程都大致过一遍才行，否则企业为什么不直接招聘一个计算机系的，而需要一个跨专业的人呢？在这种情况下，对数学系的人其实提出了很高的挑战，因为在数学系繁重的课程下，想要同时兼顾数学系和计算机系两个专业的课程是比较困难的，需要耗费巨大的时间和精力才能够做好。

balance

在这种情况下，如果是计划研究生期间转行，就需要额外花很多的功夫了。毕竟金融或者计算机其实是有一定的课程量的，不仅要保证数学系课程的学习，还要去学好金融或者计算机其实是有一定的挑战的。但是，为了最终要转行到其他领域这也是一条必经之路。因为只有补平了自己与科班人员的差距，才能够在其他领域发挥自己的数学能力。例如在计算机领域，只有编码能力过关了，才能够逐步展示数学能力。如果基础能力不太够的话，就只能够先做补基础的工作了。

experiments

数学系与其他专业的培养模式其实不太一样。在数学系，强调的一般都是听课，刷题，复习，钻研，考试等一系列流程。在这种情况下，对于想成为数学工作者的学生是比较友好的。但是对于大多数需要转行的学生而言就不算友好了。因为其他专业更多强调实践，无论是金融还是计算机都会把实践放在第一位。在金融领域，实习或者社交其实是一个比较关键的因素，但是在数学系，这些能力是无法得到锻炼的。在计算机领域，通常接触到的都是业界相对新的东西，需要在实践中一边查资料一边干活才能够逐步掌握知识点。在数学系，也会有机会进行实践，但是时间往往比其他专业会延后许多。在数学本科或者硕士阶段，其实绝大多数人是没有能力搞数学科研的，最后的毕业论文无非也就是整理一下资料和笔记而已。直到博士阶段，数学系的学生才需要进行资料的收集，论文的阅读，从而进行数学科研的工作。整体来看，数学系从事实践的时间往往会比其他专业延后许多，如果是以就业为目的来攻读数学专业其实是不合适的。

Case 3:

第三种情况，实在是不知道想学啥，如果自身条件还不错，家庭条件也还算可以的话，其实学数学专业也是一种选择。因为在数学系期间，最终的选择面还是比较广的，主要包括：

科研工作者：数学方向，金融经济方向，计算机方向等；
计算机行业；
金融行业；
教育培训行业；
其他行业。

选择面广有的时候是一件好事，有的时候不见得是一件好事。选择面太广就容易让人迷茫，尤其是对于那些目标不明确的人而言。导致的结果就是这批学生什么都去了解，什么都去学习，什么都是浅尝辄止，最终其实并没有掌握任何东西。最怕的事情就是，在数学系的时候并没有把数学弄懂，然后在其他领域又没有办法和科班竞争，从而最终荒废了自己。因此，如果实在是不清楚学什么的话，来数学系是可以的，但是肯定需要把数学弄明白。最终，突然有一天这批学生想转行或者做其他的事情，数学能力能够给人带来很大的帮助。

choice

如果萌发了转行的想法，那么最好就是找相关的人士咨询一下其他行业，或者去招聘网站上看相关的岗位需要什么样的技能，能否通过自学或者听课或者其他方式找到一份实习。然后根据自己的最终目标来反推自己当前的决定，也就是要达到这样的目标需要多少时间，多少精力才能够达成。就个人经验来看，转行是需要花大量的时间和精力的，并且越到后面（研究生或者博士阶段），转行的困难就会越来越大。因此，如果要转行的话，其实可以趁早做决定，越早决定就越有优势。

另外需要注意的一点就是，在转行期间不需要在乎数学系那一套评价指标，因为按照数学系的评价指标，不学习不科研就是不误正业了。但是对于转行的人而言，实战才是硬道理，实习才能够了解工业界的具体需求。在学校里面猛刷题或者猛看书其实只是其中的一个步骤而已，是达不到转行的要求的。

结论：

每个人都有自己的局限性，只能够根据自己的成长经历给出判断的标准和依据，并不适用于所有的人，只能建议大家根据自身的情况来选择一套适合自己的方案。毕竟人生只有一次，路也是自己走出来的。

challenge

时间序列

时间序列的联动分析

June 27, 2019 zr9558 Leave a comment

背景介绍

在互联网公司里面，通常都会监控成千上万的时间序列，用于保障整个系统或者平台的稳定性。在这种情况下，如果能够对多条时间序列之间判断其是否相关，则对于监控而言是非常有效的。基于以上的实际情况，清华大学与 Alibaba 集团在2019年一起合作了论文《CoFlux: Robustly Correlating KPIs by Fluctuations for Service Troubleshooting》，并且发表在 IWQos 2019 上。CoFlux 这个方法可以对多条时间序列来做分析，并且主要用途包括以下几点：

告警压缩和收敛；
推荐与已知告警相关的 Top N 的告警；
在已有的业务范围内（例如数据库的实例）构建异常波动传播链；

CoFlux 的整体介绍

从论文的介绍中来看，CoFlux 的输入和输出分别是：

输入：两条时间序列

输出：这两条时间序列的以下信息

波动相关性：两条时间序列是否存在波动相关性？
前后顺序：如果两条时间序列相关，那么它们的前后波动顺序是什么？是同时发生异常还是存在固定的前后顺序？
方向性：如果两条时间序列是波动相关的，那么它们的波动方向是什么？是一致还是相反？

Remark. CoFlux 的关键点就是并没有对时间序列做异常检测算法，而是直接从时间序列的历史数据（历史半个月或者一个月）出发，判断两条时间序列之间的波动相关性，并且进一步的分析先后顺序与波动方向。

从论文的介绍中来看，CoFlux 的流程图如下图所示：

coflux流程图1

如果两条时间序列 $X$ 和 $Y$ 存在波动相关性，则需要输出这两条时间序列的波动先后顺序和是否同向波动。如果两条时间序列 $X$ 和 $Y$ 并不存在波动相关性的话，则不需要判断波动先后顺序和是否同向波动。

coflux流程图2

CoFlux 的细节阐述

已知一个长度是 $n$ 的时间序列 $S=\{s_{1},\cdots,s_{n}\}$ ，对于任意一个 detector，可以得到一条关于 $S$ 的预测值曲线 $P=\{p_{1},\cdots,p_{n}\}$ 。于是针对某个 detector 可以得到一个波动特征序列 $E=\{\epsilon_{1},\cdots,\epsilon_{n}\}$ ，其中 $\epsilon_{i} = s_{i} - p_{i}$ ， $1\leq i\leq n$ 。因此，一个detector 可以对应一个波动序列特征，也是一个时间序列。因此，对于 $m$ 个 detector，可以对应 $m$ 条波动特征序列，并且它们的长度都是 $n$ 。

在 CoFlux 算法的内部，根据不同的参数使用了总共 86 个 detector，大致列举如下：

Difference：根据昨天，七天前的数据来做差分；
Holt-Winters： $\{\alpha,\beta,\gamma\} \in \{0.2,0.4,0.6,0.8\}$ ；
历史上的均值 & 历史上的中位数：1，2，3，4 周；
TSD & TSD 中位数：1，2，3，4 周；
Wavelet：1，3，5，7 天；
移动平均算法：MA，WMA，EWMA。PS：根据作者们的说法，在这里，MA等方法并不适用。

detectors

根据直觉来看，

对于任何一条时间序列 kpi，总有一个 detector 可以相对准确地提炼到其波动特征；
如果两条时间序列 $X$ 和 $Y$ 波动相关，那么 $X$ 的一个波动特征序列与 $Y$ 的一个波动特征序列应该也是相关的；

Remark. 两条时间序列的波动特征可以对齐同一个 detector，也可以不做对齐工作。如果是前者的话，时间复杂度低；后者的话，时间复杂度高。

下图是从时间序列中提取波动特征曲线的案例：

提炼时间序列的波动曲线特征只是第一步，后续 CoFlux 还有几个关键的步骤：

特征工程的扩大（amplify）：对波动序列特征进行放大，让某些波动序列特征更加明显；
Correlation Measurement：用于解决时间序列存在时间前后的漂移，两条时间序列之间存在 lag 的情况，因此需要对其中一条时间序列做平移操作；
CoFlux 考虑了历史数据（历史半个月或者一个月）作为参考，并且一个范围内的 kpi 数量不超过 60 条；

下面来一一讲解这些技术方案，对于每一条波动特征曲线（Flux-Features），按照以下几个步骤来进行操作：

Step 1：对波动特征曲线 $E=\{\epsilon_{1},\cdots,\epsilon_{n}\}$ 做 z-score 的归一化，i.e.

$\mu = \frac{\sum_{i=1}^{n}\epsilon_{i}}{n},$
$\delta = \sqrt{\frac{\sum_{i=1}^{n}(\epsilon_{i}-\mu)^{2}}{n}}.$

Step 2：对归一化之后的波动特征曲线做特征放大（feature amplification）：定义函数 $f_{\alpha,\beta}(x)$ 如下：

$f_{\alpha,\beta}(x)= \begin{cases} e^{\alpha\min(x,\beta)} - 1, \text{ when } x\geq 0,\\ -e^{\alpha\min(|x|,\beta)} + 1, \text{ when } x< 0. \end{cases}$

则 $E=\{\epsilon_{1},\cdots,\epsilon_{n}\}$ 放大之后的波动特征曲线（amplified flux feature）就是： $\hat{E}=\{f(\epsilon_{1}),\cdots,f(\epsilon_{n})\}.$

Step 3：对于两条放大之后的波动特征曲线（amplified flux features） $G=\{g_{1},\cdots,g_{\ell}\}$ 和 $H=\{h_{1},\cdots,h_{\ell}\}$ ，可以计算它们之间的相关性，先后顺序，是否同向。
令
$G_{s}= \begin{cases} \{0,\cdots,0,g_{1},\cdots, g_{\ell-s}\}, \text{ when } s\geq 0, \\ \{g_{1-s},\cdots,g_{\ell},0,\cdots,0\}, \text{ when } s< 0. \end{cases}$

这里的 0 的个数是 $|s|$ 个。其中， $-\ell<s<\ell$ 。特别地，当 $s=0$ 时， $G_{0}=\{g_{1},\cdots,g_{s}\}=G$ ，那么我们可以定义 $G_{s}$ 与 $H$ 的内积是： $R(G_{s},H) = G_{s}\cdot H,$

这里的 $\cdot$ 指的是向量之间的内积（inner product）。同时可以定义相关性（Cross Correlation）为： $CC(G_{s},H) = \frac{R(G_{s},H)}{\sqrt{R(G_{s},G_{s})\cdot R(H,H)}}.$

由于波动有可能是反向的，那么在这里我们不仅要考虑相关性是大于零的情况，也需要考虑小于零的情况。于是，

$minCC = \min_{-\ell<s<\ell}CC(G_{s},H),$
$maxCC = \max_{-\ell<s<\ell}CC(G_{s},H).$

则最小值或者最大值的指标分别是

$s_{1}=argmin_{-\ell<s<\ell}CC(G_{s},H),$
$s_{2}=argmax_{-\ell<s<\ell}CC(G_{s},H).$

令
$FCC(G,H) = \begin{cases} (minCC, s_{1}), \text{ when } |maxCC|<|minCC|, \\ (maxCC, s_{2}), \text{ when } |maxCC|\geq|minCC|. \end{cases}$

从定义中可以看出， $FCC(G,H)$ 是一个元组，里面蕴含着三个信息，分别是相关性，波动方向，前后顺序。 $FCC(G,H) \in [-1,1]$ ，越接近 1 或者 -1 就表示放大之后的波动特征曲线 $G$ 和 $H$ 越相关。正值的 $FCC(G,H)$ 表示 $G$ 与 $H$ 的波动方向相同，是正相关；负值的 $FCC(G,H)$ 表示 $G$ 与 $H$ 的波动方向想法，是负相关。通过对 $s<0$ 或者 $s\geq 0$ 的分析就可以判断先后顺序。因此，CoFlux 方法的是通过对 $FCC(G,H)$ 的分析来得到最终结果的。