目录
SLAM的数学基础（1）：什么是方差，有什么意义？
小红班上有两组同学的数学考试分数为：
第一组：小红：100分，小明：60分，小宇：20分
第二组：小蓝：70分，小华：60分，小杰：50分
那么很容易算出，第一组的平均分是60分，第二组的平均分也是60分。
这下可好，小红的100分被小宇拉了后腿。这时候，该引入一种方法，来表现这个问题。好让老师知道哪些小组的成绩差距比较大。
方差能比较好的表达一组数据离散的程度，方差大，这组数据分散的就比较大；方差小，这组数据分散的就比较小。
方差（variance）的表达公式为：
照这个公式计算，第一组的方差为：
第二组的方差为：
可以看出，第一组的方差远大于第二组。
#include <stdio.h>
float calc_variance(float samples[], int count)
{
float sum_of_samples = 0;
float average = 0;
float variance = 0;
for(int i = 0; i < count; i++)
{
sum_of_samples += samples[i];
}
average = sum_of_samples / count;
for(int i = 0; i < count; i++)
{
float temp = samples[i] - average;
variance += (temp * temp);
}
variance /= count;
return variance;
}
int main()
{
float team1[] = {100, 60, 20};
float team2[] = {70, 60, 50};
printf("variance of team 1 is %f\n", calc_variance(team1, 3));
printf("variance of team 2 is %f\n", calc_variance(team2, 3));
return 0;
}
https://www.cnblogs.com/cyberniklee/p/7923943.html
SLAM的数学基础（2）：协方差和协方差矩阵
之前我们知道，方差是一组数据的离散程度，它的公式为：
那么如果我们有几组数据，需要知道这几组数据的协同性呢？
举个例子，还是在小红，几次考试成绩如下：
入学考试：数学：80，语文：80
期中考试：数学：90，语文：70
期末考试：数学：70，语文：90
小蓝，几次考试成绩如下：
入学考试：数学：60，语文：60
期中考试：数学：70，语文：70
期末考试：数学：80，语文：80
好，我们把数据放着，先说一下概念。所谓的协方差，就是用来统计两组数据之间的协同程度，协方差矩阵是用来遍历不同组数据的方差。
协方差用公式表示为：
方差是协方差的特殊表现形式，即X=Y的时候。
OK，我们现在再来看上面的数据，根据这个公式，计算小红数学和总分的协方差为：
计算小蓝数学和总分的协方差为：
从这些数据上可以看出，小红的协方差为0，表明数学和总分没有什么关系，因为他的数学考得好，语文就考的不好。但小蓝的协方差是正的，表示他的数学进步，语文也跟着进步，总分自然提高了，是个全面发展的好青年。
那么我们想全面的表现语文、数学、总分之间的相关性，该如何表示呢？嗯，有个东西叫协方差矩阵，它可以表示多组数据之间两两之间的协方差。
比如按上面小红的数学成绩和总分的关系，形成数据集
SLAM的数学基础（3）：几种常见的概率分布的实现及验证。
分布，在计算机学科里一般是指概率分布，是概率论的基本概念之一。分布反映的是随机或某个系统中的某个变量，它的取值的范围和规律。
常见的分布有：二项分布、泊松分布、正态分布、指数分布等，下面对它们进行一一介绍。
PS：本文中谈到的PDF、PMF、CDF均为公认的缩写方式：
PDF：概率密度函数（probability density function）；
PMF：概率质量函数（probability mass function）；
CDF：累积分布函数（cumulative distribution function）。
二项分布
说起二项分布，离不开伯努利实验，二项分布就是重复N次的伯努利实验（伯努利实验，是指一种只有两种相反结果的随机试验，比如抛硬币，结果只有正面和反面；又比如投篮，只有投中和没有投中两种结果）。它的PMF可写作：
其中k为在n次实验中命中的次数，成功的概率为p。
二项分布的CDF可以写作：
例子：抛10次硬币，有2次正面朝上的概率是多少？下面分别用C++实现和用numpy证明结果
C++实现：
#include <vector>
#include <iostream>
#include <iomanip>
double calc_binomial(int n, int k, double p)
{
if(n < 0 || k < 0)
return 0.0;
std::vector< std::vector< double > > binomials((n + 1), std::vector< double >(k + 1));
binomials[0][0] = 1.0;
for(int i = 1; i < (n + 1); ++i)
binomials[i][0] = (1.0 - p) * binomials[i - 1][0];
for(int j = 1; j < (k + 1); ++j)
binomials[0][j] = 0.0;
for(int i = 1; i < (n + 1); ++i)
for (int j = 1; j < (k + 1); ++j)
binomials[i][j] = (1.0 - p) * binomials[i - 1][j] + p * binomials[i - 1][j - 1];
return binomials[n][k];
}
int main()
{
std::cout << std::fixed << std::setprecision(8) << calc_binomial(10, 2, 0.50) << std::endl;
}结果为：0.04394531
Python实现：
import numpy as np
from scipy import stats
import matplotlib.pyplot as plt
def calc_binomial():
n = 10
p = 0.5
k = 2
binomial = stats.binom.pmf(k,n,p)
print binomial
calc_binomial()结果为：0.0439453125
反之，知道投10次硬币朝上的平均概率为0.3（即平均有3次朝上），试着从10000次实验中找出规律。
用C++实现：
#include <iostream>
#include <cmath>
#include <iomanip>
#include <vector>
#include <cstdlib>
#include <ctime>
int gen_binomial_rand(int n, double p)
{
int k = 0;
for(int i = 0; i < n; i++)
{
double current_probability = ((double)rand() / (double)RAND_MAX);
if(current_probability < p)
{
k++;
}
}
return k;
}
int main()
{
srand((unsigned)time(NULL));
int gn = 10;
double gp = 0.3;
int times = 10000;
int sum_of_times = 0;
std::vector< int > result(gn);
for(int t = 0; t < times; t++)
{
int single_result = gen_binomial_rand(gn, gp);
if(single_result < gn)
{
result[single_result]++;
}
}
std::cout << std::endl;
for(int i = 0; i < gn; i++)
{
sum_of_times += result[i];
std::cout << result[i] << ",";
}
std::cout << std::endl;
std::cout << "Total: " << sum_of_times << std::endl;
return 0;
}结果为：
323,1199,2310,2631,1951,1103,367,97,18,1,
Total: 10000
拿到Python里面用图表看一下：
import numpy as np
from scipy import stats
import matplotlib.pyplot as plt
def show_binom_rvs():
n = np.array([323,1199,2310,2631,1951,1103,367,97,18,1])
plt.plot(n)
plt.show()
show_binom_rvs()显示为：
我们再来用Python的numpy和scipy的库来验证一下：
import numpy as np
from scipy import stats
import matplotlib.pyplot as plt
def calc_binom_rvs():
binom_rvs = stats.binom.rvs(n=10,p=0.3,size=10000)
plt.hist(binom_rvs, bins=10)
plt.show()
calc_binom_rvs()得到结果图片：
可以看到两次的图形包络是近似的。
泊松分布
在日常生活中，我们经常会遇到一些事情，这些事情发生的频率比较固定，但是发生的时间是不固定的，泊松分布就是用来描述单位时间内随机事件的发生概率。比如：知道一个医院平均每小时有3个小孩出生，那么下一个小时出生2个小孩的概率是多少？
泊松分布的PMF可以写作：
其中，t为连续时间长度，k为事件发生的次数，λ为发生事件的数学期望（如单位时间内发生事件的均值），e为自然底数。
就上面的例子，知道一个医院平均每小时有3个小孩出生，那么下一个小时出生2个小孩的概率是多少？用C++实现：
#include <iostream>
#include <cmath>
#include <iomanip>
double calc_poisson(int k, int lambda)
{
double result;
result = pow(lambda, k) * exp(-lambda);
int factorial = 1;
for(int i = 1; i <= k; i++)
{
factorial *= i;
}
result = result / factorial;
return result;
}
int main()
{
std::cout << std::fixed << std::setprecision(8) << calc_poisson(2, 3) << std::endl;
}结果是：0.22404181
用Python验证：
import numpy as np
from scipy import stats
import matplotlib.pyplot as plt
def calc_poisson():
lambd = 3
k = 2
y = stats.poisson.pmf(k,lambd)
print y
calc_poisson()结果是：0.224041807655
下面再用C++实现生成泊松随机数，并用Python检验：
C++实现：
#include <iostream>
#include <cmath>
#include <iomanip>
#include <vector>
#include <cstdlib>
#include <ctime>
double binary_random()
{
double rand_number = (rand() % 100);
rand_number /= 100;
return rand_number;
}
int calc_poisson(int lambda)
{
int k = 0;
double p = 1.0;
double l = exp(-lambda);
while(p >= l)
{
double r = binary_random();
p *= r;
k++;
}
return (k - 1);
}
int main()
{
int t = 10000;
int lambda = 3;
int distribution = 20;
int dist_cells[distribution] = {0};
srand((unsigned)time(NULL));
for(int i = 0; i < t; i++)
{
int n = calc_poisson(lambda);
dist_cells[n]++;
}
for(int i = 0; i < distribution; i++)
{
std::cout << dist_cells[i] << ",";
}
std::cout << std::endl;
return 0;
}运行结果为：467,1604,2298,2264,1608,952,466,217,87,29,6,2,0,0,0,0,0,0,0,0,
放入Python显示并和Python生成的比较：
import numpy as np
from scipy import stats
import matplotlib.pyplot as plt
def gen_poisson_rvs():
dist = 20
cpp_result = np.array([467,1604,2298,2264,1608,952,466,217,87,29,6,2,0,0,0,0,0,0,0,0])
py_result = np.random.poisson(lam=3,size=10000)
plt.hist(py_result,bins=dist,range=[0,dist],color='g')
plt.plot(cpp_result,color='r')
plt.show()
gen_poisson_rvs()运行并显示图表为：
指数分布
指数分布，描述的是在某一事件发生后，在连续时间间隔内继续发生的概率。比如上面的例子，知道一个医院平均每小时有3个小孩出生，刚刚已经有一个小孩出生了，那么下一个小孩在15分钟内出生的概率是多少？
指数分布的CDF可以写作：
其中t为时间长度，e为自然底数。
下面用C++实现：
#include <iostream>
#include <cmath>
#include <iomanip>
double calc_exponential(double lambda, double t)
{
double result;
result = (1 - exp((-lambda) * t));
return result;
}
int main()
{
std::cout << std::fixed << std::setprecision(8) << calc_exponential(3, 0.25) << std::endl;
}结果为：0.52763345
用Python验证：
import numpy as np
from scipy import stats
import matplotlib.pyplot as plt
def calc_expon():
lambd = 3
x = np.arange(0,1,0.25)
y = 1 - np.exp(-lambd *x)
print y
calc_expon()结果为：[ 0. 0.52763345 0.77686984 0.89460078]
其中x=0.25时结果为0.52763345
下面是生成lambda=3的指数分布随机数，样本数是10000。同样是用C++实现，Python验证。
C++实现：
#include <iostream>
#include <cmath>
#include <iomanip>
#include <vector>
#include <cstdlib>
#include <ctime>
double calc_exponential(double lambda)
{
double expon_rand = 0.0;
while(1)
{
expon_rand = ((double)rand()/(double)RAND_MAX);
if(expon_rand != 1)
{
break;
}
}
expon_rand = ((-1 / lambda) * log(1 - expon_rand));
return expon_rand;
}
int main()
{
int t = 10000;
double lambda = 3;
const int distribution = 20;
double dist_cells[distribution] = {0};
srand((unsigned)time(NULL));
for(int i = 0; i < t; i++)
{
int n = calc_exponential(lambda);
dist_cells[n]++;
}
for(int i = 0; i < distribution; i++)
{
std::cout << dist_cells[i] << ",";
}
std::cout << std::endl;
return 0;
}结果是：9515,469,16,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,
放入Python并用Python生成的对比：
import numpy as np
from scipy import stats
import matplotlib.pyplot as plt
def gen_expon_rvs():
cpp_result = np.array([9515,469,16,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0])
lambd_recip = 0.33
py_result = stats.expon.rvs(scale=lambd_recip,size=10000)
plt.plot(cpp_result,color='r')
plt.hist(py_result,color='b')
plt.xlim(0,20)
plt.show()
gen_expon_rvs()得到图片：
PS：其中最大值不一致的情形并不是计算错误，而是X轴的分布单位不一致，Python的是浮点的，而C++的代码是整形的，所以看到Python的分布最大值比较小是因为平均到了小数。
正态分布（高斯分布）
翻开任何一本讲统计的数学书，对于正态分布大抵会有相似的描述：若一个随机变量X服从一个数学期望为λ，标准差为σ的概率分布，且PDF为：
则称这个随机变量为正态随机变量。
当λ=0，σ=1时，称其为标准正态分布，PDF简化为：
在现实中，有大量的案例是符合正态分布的，比如全中国18岁以上男性人口的身高分布，170cm左右的占绝大部分，160和180的占较少部分，150以下和190以上的占极少部分。
说起正态分布的例子，有个著名的实验是不得不提的，那就是高尔顿钉板实验。
高尔顿钉板是在一块竖起的木板上钉上一排排互相平行、水平间隔相等的铁钉，并且每一排钉子数目都比上一排多一个，一排中各个钉子下好对准上面一排两上相邻铁钉的正中央。从入口处放入一个直径略小于两颗钉子间隔的小球，当小球从两钉之间的间隙下落时，由于碰到下一排铁钉，它将以相等的可能性向左或向右落下，接着小球再通过两钉的间隙，又碰到下一排铁休。如此继续下去，小球最后落入下方条状的格子内。在等可能性（即小球落在左边和落在右边的概率均为50%）的情况下，小球落下后满足正态分布。
下面的代码用C++计算模拟高尔顿实验过程，并把结果放到Python显示出来。
C++模拟实验过程：
#include <iostream>
#include <cmath>
#include <iomanip>
#include <vector>
#include <cstdlib>
#include <ctime>
int binary_random()
{
double rand_number = (rand() / (double)RAND_MAX);
if(rand_number > 0.5)
return 1;
else return 0;
}
void galton_test(int num_of_cells, int num_of_balls)
{
srand((unsigned)time(NULL));
std::vector< int > cells(num_of_cells);
int rand;
for(int i = 1; i <= num_of_balls; i++)
{
int cell = 0;
for(int j = 1; j < num_of_cells; j++)
{
int rand = binary_random();
cell += rand;
}
cells[cell]++;
}
std::cout << std::endl;
for(int i = 0; i < num_of_cells; i++)
{
std::cout << cells[i] << ",";
}
std::cout << std::endl;
}
int main()
{
galton_test(20, 5000);
}结果为：0,0,3,14,45,95,264,481,720,886,911,741,452,240,95,45,6,1,1,0,
结果每次都不一样，但将其放入以下Python代码显示：
import numpy as np
from scipy import stats
import matplotlib.pyplot as plt
def show_galton():
n = np.array([0,0,3,14,45,95,264,481,720,886,911,741,452,240,95,45,6,1,1,0])
plt.plot(n)
plt.show()
show_galton()得到以下图片：
可以看出，是个明显的钟型曲线，说明高尔顿实验是满足正态分布的。在这个实验中，我们还可以去调整num_of_cells, num_of_balls的值，可以看出，当num_of_cells的值越大，曲线越陡峭，越小，曲线越平坦；num_of_balls的值越大，曲线就越像正态分布。说到这里可能大家就会想的到，这个实验中，num_of_cells的值可以认为就是正态分布的σ，而我们设定的随机数0，1会影响到正态分布的λ，有兴趣的朋友可以改一下上面的例程，将随机数生成改为大于0.5，或小于0.5，可以观察到最后的曲线中轴线会向左偏和向右偏。
说到这里，下面的公式应该是理所当然的了：
若一个随机变量X服从一个数学期望为λ，尺度参数为σ的概率分布，记作
说到这里，通过上面的实验大家应该大概知道高斯分布是个什么东西了。那么如何编程实现生成符合高斯分布的随机数呢？
生成高斯分布的随机数有多种方法：
（1） Box-Muller变换算法
（2） 利用中心极限定理迭代法
（3） Ziggurat算法
等。其中，在效率和通用性方面比较均衡的是Box-Muller算法，C++11和Python的数学库里面基本都是用的它。
它的原理是：
随机抽出从[0,1]中符合均匀分布的数a和b，然后令：
那么这两个数都是符合正态分布的。
若想产生服从期望是λ，标准差是σ的正态分布，那么：
下面，我们用C++来实现：
#include <iostream>
#include <cmath>
#include <iomanip>
#include <vector>
#include <cstdlib>
#include <ctime>
#include <limits>
double calc_gaussian(double sigma, double lambda)
{
static const double epsilon = std::numeric_limits<double>::min();
static const double two_pi = (2.0 * 3.14159265358979323846);
static double z1;
static bool generate;
generate = !generate;
if (!generate)
return z1 * sigma + lambda;
double a, b;
do
{
a = rand() * (1.0 / RAND_MAX);
b = rand() * (1.0 / RAND_MAX);
}while ( a <= epsilon );
double z0;
z0 = sqrt(-2.0 * log(a)) * cos(two_pi * b);
z1 = sqrt(-2.0 * log(a)) * sin(two_pi * b);
return z0 * sigma + lambda;
}
int main()
{
int t = 10000;
double lambda = 10;
double sigma = 1;
const int distribution = 20;
double dist_cells[distribution] = {0};
srand((unsigned)time(NULL));
for(int i = 0; i < t; i++)
{
int n = calc_gaussian(sigma, lambda);
dist_cells[n]++;
}
for(int i = 0; i < distribution; i++)
{
std::cout << dist_cells[i] << ",";
}
std::cout << std::endl;
return 0;
}得到结果：0,0,0,0,0,0,12,190,1374,3372,3519,1325,196,12,0,0,0,0,0,0,
放入Python显示：
SLAM的数学基础（4）：先验概率、后验概率、贝叶斯准则
贝叶斯公式
假设有事件A和事件B，可以同时发生但不是完全同时发生，如以下韦恩图所示：
其中，A∩B表示A和B的并集，即A和B同时发生的概率。
如此，我们很容易得出，在事件B发生的情况下，事件A发生的概率为：
这个P(A|B)就是条件概率(Conditional Probability)。
同理，在事件A发生的情况下，事件B发生的概率为：
由以上式子可得：
再调整一下，变成：
这个就是著名的贝叶斯公式的基本形态了，其中：
P(A|B)叫做后验概率(Posterior Probability)
P(A)叫做先验概率(Prior Probability)
P(B|A)/P(B)叫做似然度(Likelihood)
那么我们可以看出，贝叶斯定理可以比较简单的归纳为：
后验概率=先验概率*似然度
在日常使用中，如贝叶斯分类、贝叶斯回归、贝叶斯滤波等算法，普遍使用迭代和归一化的方法来计算似然度
全概率公式
在日常使用中，如贝叶斯分类、贝叶斯回归、贝叶斯滤波等算法，普遍使用迭代和归一化的方法来计算似然度，为了更好的了解归一化的方法，这里还有一个基础概念，叫全概率公式。
接着之前的说法，很明显有
OK，先举几个经典的例子来帮助理解：
一座别墅在过去的 20 年里一共发生过 2 次被盗，别墅的主人有一条狗，狗平均每周晚上叫 3 次，假设在盗贼入侵时狗叫的概率被估计为 0.9，问题是：在狗叫的时候发生入侵的概率是多少？
我们假设：
事件A：狗在晚上叫
事件B：盗贼入侵
以天为单位统计,有：
可以看到，最终的结果基本是不可能。
这只是一个最简单的例子，其中，别墅在过去的20年中被盗2次、狗平均每周晚上叫3次这些都是先验统计，而在现实应用中往往是根据上一步的事件去推断下一步的事件，这个迭代的过程往往是很多算法实现的基础。
下面，我们举一个更复杂一些的例子：
假设有两个各装了100个球的箱子，A箱子中有70个红球，30个蓝球，B箱子中有30个红球，70个蓝球。假设随机选择其中一个箱子，从中拿出一个球记下球色再放回原箱子，如此重复12次，记录得到8次红球，4次蓝球。问题来了，你认为被选择的箱子是A箱子的概率有多大？
我们假设得到小球的结果顺序如下：红红红红红红红红蓝蓝蓝蓝，用C++实现：
#include <iostream>
#include <cmath>
#include <iomanip>
#include <vector>
#include <cstdlib>
#include <ctime>
#include <limits>
class balls_case_problem
{
private:
double p_a;
double p_b;
double p_red_in_a;
double p_blue_in_a;
double p_red_in_b;
double p_blue_in_b;
public:
balls_case_problem(int red_balls_in_a, int blue_balls_in_a, int red_balls_in_b, int blue_balls_in_b)
{
p_a = 0.5;
p_b = 1 - p_a;
p_red_in_a = (double)red_balls_in_a / (double)(red_balls_in_a + blue_balls_in_a);
p_blue_in_a = 1 - p_red_in_a;
p_red_in_b = (double)red_balls_in_b / (double)(red_balls_in_b + blue_balls_in_b);
p_blue_in_b = 1 - p_red_in_b;
}
~balls_case_problem(){};
void got_red()
{
p_a = (p_red_in_a * p_a) / ((p_red_in_a * p_a) + (p_red_in_b * p_b));
p_b = 1 - p_a;
}
void got_blue()
{
p_a = (p_blue_in_a * p_a) / ((p_blue_in_a * p_a) + (p_blue_in_b * p_b));
p_b = 1 - p_a;
}
double get_p_a()
{
return p_a;
}
double get_p_b()
{
return p_b;
}
};
int main()
{
int red_ball_results[] = {1, 1, 1, 1, 1, 1, 1, 1, 0, 0, 0, 0};
balls_case_problem problem(70, 30, 30, 70);
for(int i = 0; i < sizeof(red_ball_results) / sizeof(int); i++)
{
if(red_ball_results[i] == 1)
{
problem.got_red();
}else{
problem.got_blue();
}
std::cout << "Probility of chose case A: " << problem.get_p_a() << std::endl;
}
return 0;
}运行结果为：
Probility of chose case A: 0.7
Probility of chose case A: 0.844828
Probility of chose case A: 0.927027
Probility of chose case A: 0.967365
Probility of chose case A: 0.985748
Probility of chose case A: 0.993842
Probility of chose case A: 0.997351
Probility of chose case A: 0.998863
Probility of chose case A: 0.997351
Probility of chose case A: 0.993842
Probility of chose case A: 0.985748
Probility of chose case A: 0.967365
被折叠的 条评论
为什么被折叠?
到【灌水乐园】发言
