中科大-算法分析与设计-概率算法复习知识点
二、概率算法
概率算法分类：
- 数字算法：求近似解
- Monte Carlo 算法：解判定问题或只有一个正确解的问题。答案未必正确。
- Las Vegas 算法：决不返回错误答案。
- Sherwood 算法：总是给出正确答案。
2.1 数字算法
找到数字问题的近似解
\(\Pi\)值计算：投飞镖模拟，单位圆和外切正方形。
\[\frac{落在圆内}{总数} = \frac{k}{n} =\frac{\pi r^2}{4r^2} = \frac{\pi}{4} \\ \Rightarrow \pi = 4 \frac{k}{n} \]Darts(n){ k=0; for i=1 to n do { x = uniform(0,1); y = uniform(0,1); if(x*x+y*y<=1) then k++; } return 4k/n; }
数字积分：Monte Carlo 积分，但不是本课中的 MC 算法。
计算积分\(S=\int_0^!f(x)dx\)
HitorMiss 算法：单位面积正方形投飞镖，落入积分函数下方的为 k，\(\frac{k}{n}=\frac{S}{1}\)
HitorMiss(f, n){ k=0; for i=1 to n do { x = uniform(0,1); y = uniform(0,1); if(y<=f(x)) then k++; } return k/n; }Crude 算法：积分区间上随机取点，取算术平均值乘区间宽度。\(\int_a^bf(x)dx=(b-1)\frac{1}{n}\sum_{i=1}^{n}f(x_i)\)
Crude(f,n,a,b){ sum = 0; for i=1 to n do{ x = uniform(a,b); sum = sum + f(x); } return (b-a)*sum/n; }相对而言，Crude 算法方差更小，但是 HitorMiss 算法在相同时间迭代次数更多。
解题思路：\(\pi=4\int_0^1\sqrt{1-x^2}dx\)
Trapezoid 确定性算法：梯形算法，将区间划分为 n-1 个子区间
\(S=\delta*(f(a+\delta)+f(a+2\delta)+\cdots+\frac{f(a)+f(b)}{2})\)
Trapezoid(f,n,a,b){ // 假设n>=2 delta = (b-a)/(n-1); sum = (f(a)+f(b))/2; for(x=a+delta; x<=b-delta; x+=delta){ sum = sum + f(x); } return sum * delta; }确定性算法有时求不出解，概率算法适用于高维定积分。
概率计数：估计整数值，如集合的大小或者集合中不同对象数目的估计。
求集合的大小：有放回的随机抽取元素，k 是第 1 次出现重复元素之前所选出的数目。
\(n足够大时,k=\beta\sqrt{n}=\sqrt{\frac{\pi}{2}n},即n=\frac{2k^2}{\pi}\)
SetCount(X){ k = 0; S = {}; a = uniform(X); do{ k++; S = S U {a}; a = uniform(X); }while(a不属于S) return 2k**2/pi }多重集合不同对象数目估计：将单词映射到长度为 m 的位串，\(\pi(y,b)=\begin{cases}i,y[i]=b的最小i \\ k+1,y[i]!=b \end{cases}\)
WordCount(){ y[1...m+1] = 0; for 磁带上的每个单词 x do { i = pi(h(x), 1); // x的散列值中等于1的最小位置，表示x是以 00...01....开头的 y[i] = 1; } return pi(y, 0); // 返回y中等于0的最小位置 } 2^(k-2)<=n<=2^k，更准确的n=2^k/1.54703
2.2 Sherwood 算法
总是给出正确的答案，平滑不同输入实例的执行时间。
用途：平均性能较优，但是最坏性能较差时才使用。比如选择和排序，对元素的相对顺序敏感。
基本过程：
- 将被解实例变换到一个随机实例。（预处理）
- 用确定算法解随机实例，得到一个解。
- 将此解变换为原实例的解。（后处理）
RH(x){ // 用Sherwood算法计算f(x) n = length(x); // x的size为n r = uniform(A[n]); // 随机选取一个元素 y = u(x,r); // 将原实例x转化为随机实例y s = f(y); // 用确定算法求y的解s return v(r,s); // 将s的解还原为x的解 }
选择和排序的 Sherwood 算法：将输入实例打乱，无需后处理。在调用确定性算法前，先洗牌
Shuffle(T){
n = length(T);
for i=1 to n-1 do {
// 在T[1...n]中随机选取一个元素放在T[i]上
j = uniform(1...n);
sqap(T[i],T[j]);
}
}
离散对数计算：\(a=g^x mod \ p, 计算x=log_{g,p}a\)
简单算法：
- 计算\(g^x mod \ p,0 \le x\le p-1\)，x 有限个取值，因为离散对数\(<g>\)是个循环群
- 验证\(a=g^x mod \ p\)是否成立
dlog(g,a,p){ // 当这样的算法不存在时，算法返回p x = 0; y = 1; do{ x++; y = y*g; // 计算y=g**x }while(a!=y%p && x!=p); return x; }定理：
- \(log_{g,p}(st \% p)=(log_{g,p}s+log_{g,p}t)\%(p-1)\)
- \(log_{g,p}(g^r \% p)=r,0\le r \le p-2\)
Sherwood 改造算法：
dlogRH(g, a, p){ r = uniform(0..p-2); b = g^r mod p; // 求幂模b=g**r%p c = ba mod p; // 定理1 y = dlog(g,c,p); // 使用确定性算法求log_{g,p}c, y=r+x return (y-r) mod (p-1); // 求x }\[\begin{align} y &= log_{g,p}(ba \% p) = (log_{g,p}b+log_{g,p}a)\%(p-1) \\ &= (log_{g,p}(g^r \% p)+log_{g,p}a)\%(p-1) \\ &= (r+x)\%(p-1) \end{align} \]
搜索有序表：数值数组和指针数组构成有序的静态链表。
时间为\(O(n)\)的确定算法：直接顺序搜索，平均时间\(O(\frac{n}{2})\)
Search(x, i){ while x > val[i] do i = ptr[i]; return i; } A(x){ return Search(x, head); }时间为\(O(n)\)的概率算法：从随机的一个点左右开始顺序搜索，平均时间\(O(\frac{n}{3})\)
D(x){ i = uniform(1..n); y = val[i]; case{ x<y: return Search(x, head); // 从头开始搜索 x>y: return Search(x, ptr[i]); // 从下一个位置开始搜索 otherwise: return i; // case3 x==y } }平均时间为\(O(\sqrt{n})\)的确定算法：在前\(\left \lfloor \sqrt{n} \right \rfloor\)个数中找不大于 x 的最大整数的下标，从该数开始向后查找。
B(x){ // 设x在val[1..n]中 i = head; max = val[i]; // max初值是表val中的最小值 for j=1 to [sqrt(n)] do { // 在前根号n个数中找不大于x的最大整数y的下标 y = val[j]; if max < y <= x then { i = y; max = y; } } return Search(x, i); }
2.3 Las Vegas 算法
特点：
- 更有效率的算法
- 时间上界可能不存在
- 可能找不到解，但如果找到一定是正确的
- 成功的概率随执行时间的增加而增加
算法的一般形式：\(LV(x,y,success),x为输入,y为输出\)
约定记号：
- \(p(x),算法成功的概率\)
- \(s(x),算法成功时的期望时间\)
- \(e(x),算法失败时的期望时间\)
一般会多次运行 LV 算法，知道成功
Obstinate(x){ repeat LV(x,y,success); until success; return y; }找到正确解的期望时间：
\[\begin{align} t(x) &= p(x)s(x)+(1-p(x))(e(x)+t(x)) \\ &=s(x)+\frac{1-p(x)}{p(x)}e(x) \end{align} \]
2.3.1 八皇后问题
要求：
- 行不冲突：行号相等
- 列不冲突：列号相等
- 主对角线不冲突：行列号之差相等
- 辅对角线不冲突：行列号之和相等
传统回溯法：分别决定每一行放在哪个位置
i = j = 1;
while i ≤ 8 do { // 当前行号i≤ 8
检查当前行i：从当前列j起向后逐列试探，寻找安全列号；
if 找到安全列号 then {
放置皇后，将列号记入栈，并将下一行置为当前行(i++);
j=1; //当前列置为1
} else {
退栈回溯到上一行，即i--;
移去该行已放置的皇后，以该皇后所在列的下一列作为当前列;
}
}
Las Vegas 算法：
约定：
- \(try[1...8]表示第i个皇后放在(i,try[i])上\)
- \(try[1...k]称为k-promising\)，当且仅当满足 k 皇后要求
- \(\forall i!=j \in [1,k],try[i]-try[j]\notin \{i-j,0,j-i\}\)
基本步骤：
- 遍历改行所有的列位置是否可行，\(\frac{1}{nb}\)概率下选用该列
- 如果有解则放，没解则退出
QueensLV(success){
// 贪心的LV算法，所有皇后都是随机放置
col, diag45, diag135 = {}; // 列、两个对角线初始化为空集
k = 0; // 行号
repeat
nb = 0; // 计数器，nb值为(k+1)th皇后的open位置总数
// 遍历列号，试探(k+1,i)是否安全
for i=1 to 8 do {
if (i!∈col) and (i-k-1!∈diag45) and (i+k+1!∈diag135) then {
// 列i对(i+1)th皇后可用，但不马上放置
nb = nb + 1;
if uniform(1..nb)=1 then // 或许放在第i列
j = i;
}
}
if(nb>0) then {
// nb=0时无安全位置
// 在所有nb个安全位置上，(k+1)th皇后选择位置j列的概率为1/nb
k = k+1;
try[k] = j;
col = col ∪ {j};
diag45 = diag45 ∪ {j-k};
diag135 = diag135 ∪ {j+k};
}
until nb=0 or k=8; // 当找不到合适的位置，或者try是8-promising时，退出
success = (nb>0);
}
算法 12 行，在所有可能的列位置上，选中某一列的概率都是\(\frac{1}{nb}\)，从最后一项往前倒推。
改进算法：联合回溯法和LV 算法，指定使用 LV 算法只确定前几行的位置。
...上述算法...
until nb=0 or k=stepVegas; // 已随机放好stepVegas个皇后
if nb>0 then
backtrace(k, col, diag45, diag135, success);
else
success = false;
2.3.2 模 p 平方根
\(x \equiv y^2 (mod \ p),x\in[1,p-1],p为奇素数\)，则
- x 为模 p 的二次剩余
- y 为 x 模 p 的平方根
定理：
- 任何一个模 p 的二次剩余有两个不同的平方根 \(b,(p-b)\)
- 1 到 p-1 之间的整数恰有一半时模 p 的二次剩余
- \(\forall x \in [1,p-1],p为奇素数,有x^{(p-1)/2} \equiv \pm1(mod \ p)\)
- \(x是模p的二次剩余，当且仅当x^{(p-1)/2} \equiv 1(mod \ p)\)
如何判定x 是否为模 p 的二次剩余？计算\(x^{\frac{p-1}{2}}(mod \ p)\)是否为 1。
如何计算x 模 p 的两个平方根？
-
\(if \ p\equiv 3 (mod \ 4)\)，则平方根为 \(\pm x^{(p+1)/4}\)
-
\(if \ p\equiv 1 (mod \ 4)\)，只能使用 LV 算法
例子：\(p=53\equiv1(mod \ 4),x=7\)，求 7 模 53 的平方根。
基本步骤：
- 随机取\(a\in[1,p-1]\)
- 计算 c 和 d，\((a+\sqrt{x})^{(p-1)/2} \equiv c+d\sqrt{x} \ (mod \ p)\)
- 若\(d=0\)，则返回假
- 若\(c=0\)，则找到 y，\(dy \equiv 1 \ (mod \ p)\)
rootLV(x, p, y, success){ // 计算y
a = uniform(1..p-1); // 并不知道a取多少
if a**2 % p == x % p then { // 可能性很小
success = true; y =a;
}else{
计算c和d，使得0<=c,d<=p-1，(a+sqrt(x))**((p-1)/2) % p == c+d*sqrt(x) % p;
if d==0 then success = false;
else{ // c=0
success = true;
计算y使得，d*y % p == 1, 1<=y<=p-1 // 修改Euclid算法可求
}
}
}
算法成功的概率大于 0.5，因此平均调用两次即可求得 x 的平方根。
2.3.3 整数的因数分解
整数的素因数分解，\(n=p_1^{m_1}p_2^{m_2}\cdots p_k^{m_k}\)
朴素 split 算法：找到 n 的一个非平凡因数
split(n){
// n是素数返回1，否则返回找到的n的最小非平凡因数
for i=2 to [sqrt(n)] do
if n % i == 0 then return i;
return 1;
}
合数的二次剩余：上节的二次剩余要求 p 为奇素数，则恰有两个不同的平方根。若\(n=pq,p和q为不同素数\)，则二次剩余恰有4 个平方根。
Dixon 因数分解算法：
基本步骤：
- 找到与 n 互素的整数 a 和 b，使得 \(a^2 \equiv b^2 \ (mod \ n) 且 a !\equiv \pm b \ (mod \ n)\)
- 则 n 和 a+b 的最大公因子是一个非平凡因子，\(x=gcd(n, a+b)\)。
Dixon(n,x,success){ // 找合数n的某一非平凡因子
if n是偶数 then {
x = 2; success = true;
}else{
for i=2 to [log_{2}n] do {
if n**(1/i) 是整数 then {
x = n**(1/i); success = true; return;
}
}
找到a, b，使得 a**2 % n == b**2 % n;
if a % n == +-b % n then success = false;
else {
x = gcd(a+b, n); // 最大公因子
success = true;
}
}
}
k-平滑：一个整数 x 的所有素因子均在前 k 个素数中。
确定 a 和 b 的取值：提前选定 n 和 k。
随机选择\(x\in[1,n-1]\)。
- 若 x 不与 n 互素，则找到了非平凡因子
- 否则，\(y=x^2 \ mod \ n\)，若 y 是 k 平滑的，则将 x 和 y 的因数分解保存在表里。
- 一直重复，直到选择了k+1 个互不相同的整数 x。
在 k+1 个因式分解式中找到一个非空子集\(\{y_1y_2\cdots y_k\}\)，使得相应的因式分解的积中前 k 个素数的指数均为偶数（包含 0），如\(y_1y_2y_4y_8=2^63^45^47^011^213^217^2\)
令\(S=\{y_{s_1},y_{s_2},\cdots,y_{s_k}\},a=\prod_{y_i\in S} x_i,b=\sqrt{\prod_{y_i\in S} y_i}\)，\(if \ a !\equiv \pm b \ (mod \ n)\)，则只需要求 a+b 和 n 的最大公因子即可。
时间分析：k 的取值将会影响 \(x^2 \ mod \ n\)是 k 平滑的可能性，以及因数分解 y 时的时间。通常取如下时间
2.4 Monte Carlo 算法
存在某些问题无法找到有效的算法获得正确解。MC 以高概率找到正确解。
约定：
- p-正确：一个 MC 算法以不小于 p 的概率返回一个正确的解
- 相容的、一致的：一个 MC 算法对同一实例返回相同的正确解。
多次调用同一个算法，选择出现次数最多的解，以增加一致的，p-正确算法成功的概率。
有偏算法：
- 偏真算法：MC 算法返回 true 时的解一定正确，只有返回 false 时才有可能产生错误的解。
- 偏\(y_0\)算法：算法返回\(y_0\)时总是正确的，返回非\(y_0\)时 p 概率正确。
重复调用一个一致的，p-正确的，偏\(y_0\)的 MC 算法 k 次，得到一个 \(1-(1-p)^k\) -正确的算法。
2.4.1 主元素问题
一个具有 n 个元素的数组，若某一元素的个数大于\(\frac{n}{2}\)，则其为主元素。
判断 T 是否有主元素：随机设置一个 x，然后看看有多少个
maj(T){
i = uniform(1..n);
x = T[i];
k = 0;
for j=1 to n do
if T[j] == x then
k = k+1;
return k>n/2;
}
偏真，0.5-正确的算法
重复调用降低错误率：
maj2(T){
if maj(T) then return true; // 1次成功
else return maj(T); // 调用2次
}
偏真，0.75-正确的算法
算法改进：控制算法出错概率小于\(\varepsilon > 0\)，应调用次数为\(k=\left \lceil lg\frac{1}{\varepsilon} \right \rceil\)
// e 为错误概率
majMC(T, e){
k = [lg(1/e)]; // 上取整
for i=1 to k do
if maj(T) then return true;
return false;
}
// O(nlg(1/e))
2.4.2 素数测定
判定一个整数是否为素数，尚未找到有效的确定性算法或 LV 算法
简单的概率算法：
prime(n){
d = uniform(2..[sqrt(n)]); // 下取整
return (n mod d)!=0;
}
Fermat 小定理：\(if \ n 是素数,则\forall a \in [1,n-1],有a^{n-1} mod \ n = 1\)
变换命题：\(if \ \exist a \in [1,n-1],使a^{n-1}mod \ n != 1, 则n不是素数\)
素性测定算法：使用 Fermat 小定理，偏假的。即算法返回假，一定为假。
Fermat(n){
a = uniform(1..n-1);
if a**(n-1) mod n == 1 then return true; // 未必正确
else return false; // 一定正确
}
伪素数和素性伪证据：符合 Fermat 小定理的合数，称为以 a 为底的伪素数，a 称为 n 的素性伪证据。
将 Fermat 测试重复多次，也不能降低误差。
强伪素数：n 是大于 4 的奇数，s 和 t 是整数，使得\(n-1=2^st,t为奇数\)，集合\(B(n)\)满足如下条件：
- n 为素数，则 \(\forall a \in [2,n-2], a \in B(n)\)
- n 为合数，\(if \ a \in B(n)\)，则称 n 为一个以 a 为底的强伪素数，a 为 n 素性的强伪证据。
Btest(a, n){
// n为奇数，2<=a<=n-2, 返回true则代表a属于B(n)。则n为强伪素数或素数
s=0; t=n-1; // t开始为偶数
repeat
s++;t/=2;
until t mod 2 = 1; //n-1=(2**s)t, t为奇数
x = a**t mod n;
if x=1 or x=n-1 then return true; // 满足条件1或2
for i=1 to s-1 do{
x = x**2 mod n;
if x=n-1 then return true; // 满足条件2
}
return false;
}
0.75-正确
强伪证据的数目比伪证据数目少很多。
\(强伪证据 \Rightarrow 伪证据\)
Miller-Rabin 测试：0.75-正确，偏假的
MillRab(n){ // 奇n>4，返回真时表示素数，假时表示合数
a = uniform(2..n-2);
return Btest(a,n); // 测试n是否为强伪素数
}
RepeatMillRob：重复实验 k 次，\(1-4^{-k}\)-正确的 MC 算法
RepeatMillRob(n,k){
for i=1 to k do
if MillRob(n) == false then return false; // 一定是合数
return true;
}
若给出出错概率不超过\(\varepsilon\)，则重复次数为 \(k=\left \lceil \frac{1}{2}lg\frac{1}{\varepsilon} \right \rceil\)
确定 n 素性时间为：\(O(lg^3n \ lg\frac{1}{\varepsilon})\)
2.4.3 矩阵乘法验证
验证\(矩阵AB=C\)
设 \(X_{1 \times n}=(a_1,a_2,\cdots,a_n),a_i\in\{0,1\}\)，则问题转化为 \(XAB=XC\)
基本步骤：
- 计算\(XA\)
- 计算\(XA\)与\(B\)的乘积
- 计算\(XC\)
GoodProduct(A,B,C,n){
for i=1 to n do
x[i] = uniform(0..1);
if XAB=XC then return true;
else return false;
}
偏假的，0.5-正确
返回假一定为假
改进算法：重复 k 次
RepeatGoodProduct(A,B,C,n,k){
for i=1 to k do // 重复k次
if GoodProduct(A,B,C,n) == false then return false;
return true;
}
偏假的，\(1-2^{-k}\)-正确
给出错误概率：即$2^{-k}=\varepsilon $
GP(A,B,C,n,e){
k = [lg(1/e)]; // 上取整
return RepeatGoodProduct(A,B,C,n,k);
}
// O(n**2log(1/e))
本文来自博客园，作者：小陆斑比，转载请注明原文链接：https://www.cnblogs.com/bamlubi/p/15763708.html
