抽样方法(2)--拒绝抽样和高斯分布

在上篇博客抽样方法(1)--随机数和简单随机抽样中,我介绍了如何生成随机数和均匀分布。进而通过均匀分布我们可以用Inverse Transform Method对一些简单分布进行抽样,但在博客的最后,我提到到由于高斯分布的累积分布函数是不可积的,因此用 Inverse Transform Method 不能直接对高斯分布进行抽样。这篇博客,我们就来为高斯分布寻找一个抽样方法。在继续后面的内容之前,注意到由于$N(\mu,\delta^2)=\mu+\delta N(0,1)$,因此要生成均值为$\mu$,方差为$\delta^2$的高斯分布,只需要我们能生成标准的高斯分布即可。针对前面高斯分布抽样的存在问题,我们可以从两个方向入手:

  1. 既然高斯分布的CDF: $F(x) = \int_{-\infty}^x \frac{1}{\sqrt{2\pi}} exp(-\frac{1}{2} x^2) dx$是不可积的,那么我们可以通过一些手段使得上式可积或者得到这个积分的一个近似。

  2. 采用其他不需要计算累计分布函数的抽样方法。

阅读更多…

C++ 对象模型

我觉得学习中最幸福的时刻是 "aha moment": 当你弄懂一个思考很久的知识点,然后发现以前郁结于心的很多问题都解释得通了的那一刻。在我学习C++的经历中,印象最深的"aha moment"当属弄懂对象模型的时刻,因此当我想写关于C++的博客时,我决定利用"贪婪算法":先写我觉得最有趣的东西。

本篇博客中的例子都是用 clang3.5在 MacOSX 10.11 64位下编译得到,有些例子在不同的系统和编译器下可能结果不同

C++ 中最简单的对象类型是 Plain Old Data (POD)类型,例如下面这个结构体:

struct toy{
    int a;
    int b;
}

这样的结构体和C语言的结构体别无二致,它没有构造函数和析构函数,它包含的数据都是原始类型,它也没有继承别的类。像上面的这种C++结构体是和C兼容的,可以用C语言的内存函数如memcpy,memset来操作这个结构体。它的内存布局也相对简单,结构体的大小就是结构体内元素的大小加上内存对齐。比如上面的toy这个结构体,用sizeof(toy)可以得到它的大小为8(Bytes),即两个int型变量的大小的和,如果把toy里b的类型改为char型,虽然一个int加上一个char大小是5,但是由于内存对齐的要求,仍然要为toy分配8字节的内存,C++标准要求类中的数据布局和变量的声明顺序一致,在toy 中,由于a先声明,b后声明,在toy的实例中,a也要在前,b在后,地址由低到高。这里面隐藏的一个问题是,由于a占了前4个字节,b 只占一个字节,它在后面4个字节任何位置都是对齐的,那么编译器会把b放在哪里呢?可以用下面的程序验证:

#include<iostream>
struct toy{
    int a;
    char b;
};
int main (){
    toy t;
    t.b = 'a';
    char* p = reinterpret_cast<char*>(&t);
    for(size_t i=0;i<4;++i)
        std::cout<<i+5<<":"<<p[i+4]<<" ";
    return 0;

 /*输出:
    5:a 6: 7: 8:
 */
}

上面的程序输出了t从第5个字节起的所有值,注意 6,7,8内存位置都是不可打印的值,5号位置是程序main函数第二行写入b的值。也就是说b的内存是在后四个字节的起始位置。另一种情况是,如果在toy中再增加一个int16_t 类型的数据,它是在第6,7个字节吗?注意这个时候由于它是需要对齐的,它的大小是2个字节,因此它应该对齐在第7个字节的内存位置。第6号位置就被跳过了,这时的toy内存布局为:

address:| 0 | 1 | 2 | 3 |  4  | 5 | 6 | 7 |
    toy{|      int      | char|   |int16_t|}

阅读更多…

抽样方法(1)--随机数和简单随机抽样

《统计模拟》是 M.ROSS 写的一本小册子,由于这本小册子携带很方便,我通常是在零碎时间看,这样断断续续看了一个学期,导致我学的相当不系统。于是我决定在博客里记录一下书中的内容,帮助自己思路的同时也回顾一下机器学习中的各种 sampling 方法。想想内容其实还是挺多的,所以我初步计划分成三个部分写。这篇博客其中的第一部分,主要内容是(伪)随机数的生成算法,以及对一些简单的分布进行抽样的方法。

1. 随机数的生成算法

在介绍具体算法前,先解释一下,我们在各种编程语言中用类似random()这样的函数调用得到的随机数大多数是伪随机数,它们通常都是给的一个初始值(称为随机数种子),然后按照一个给定的公式计算下一个值。线性同余发生器是其中最简单的之一,它的计算公式是 $ X^{n+1} = AX^{n} + B \mod M $ 公式中的A,B,M都是事先给定的常数,而$X_0$就是我们的随机数种子。从这个公式能得到的一个很直观的结论就是这个算法最多产生$N (N\leq M)$个随机数后就会产生重复的随机数序列,这个就N叫做这个算法的周期。通常我们希望这个算法的周期越长越好,即N越大越好,但是在这个算法里,N要达到$M$也是不容易的,它对$A,B,M$的选取有特殊的规定,比如要求$B$和$M$互质以及一些其他条件。虽然这个算法简单,但是它的优点很明显,就是计算量小,每生成一个随机数只需要一次乘法一次加法一次取模运算。因此这个算法通常用在一些追求速度的算法中,比如 Word2Vec 中生成随机数就是用的线性同余发生器,它的表达式是 next_random = next_random*25214903917 +11,可以看到它的A和B选择了两个很fancy的质数,这里的 next_random C语言中的 long long 型的变量,所以这里的M就隐含着是long long 型变量的最大值。这个算法参数的随机数是否均匀分布可以用下面这段代码来验证,这段代码从这个算法中生成了N个随机数,然后画出了这些随机数归一化后的累积概率分布和概率密度函数。