投稿
  1. 高性能架构探索首页
  2. C/C++

概率与期望

admin • • C/C++ • 阅读 176

定义

  • 概率，就是某个随机事件出现的可能性大小。

  • 若 \(X\) 是一个离散型的随机变量，可能值为 \(x_1,x_2…\)，对应的概率分别为 \(p_1,p_2…\)，那么它的期望值为 \(E(x)=\sum_i \limits p_ix_i\)。

期望的线性性

\[E(x+y)=E(x)+E(y)
\]

证明：

\[E(x+y)=\sum_i \sum_j(i+j)*P(i=x,j=y)
\]

\[=\sum_i\sum_ji*P(i=x,j=y)+\sum_i\sum_jj*P(i=x,j=y)
\]

\[=\sum_ii*P(i=x)+\sum_jj*P(j=y)
\]

\[=E(x)+E(y)
\]

进而有：

\[E(ax)=aE(x)\\
\]

\[E(\sum_i a_i x_i)=\sum_i a_i E(x_i)\\
\]

当随机变量 \(x\) 和 \(y\) 独立时，有：

\[E(xy)=E(x)E(y)
\]

条件概率

我们记 \(p(A|B)\) 表示在已知事件 \(B\) 发生时事件 \(A\) 发生的概率，条件概率可以用以下公式计算：

\[P(A|B)=\frac{P(AB)}{P(B)}
\]

其中 \(p(AB)\) 表示事件 \(B\) 和事件 \(A\) 同时发生的概率，\(p(B)\) 表示事件 \(B\) 发生的概率。

贝叶斯公式

由条件概率的计算方法，我们容易得到贝叶斯公式:

\[P(A|B)=\frac{P(B|A)P(A)}{P(B)}
\]

全概率公式

如果随机变量 \(x\) 有 \(k\) 个取值，分别为 \(x_1, x_2,\ldots ,x_k\)，那么：

\[p(A)=\sum^{k}_{i=1} {p(A|x=x_i)p(x=x_i)}
\]

[CTSC2017]游戏

根据期望的线性性，我们可以单独计算每个未知的局面的获胜概率。

设 \(x\) 是一个未知局面，显然 \(x\) 获胜的概率只与 \(x\) 左边的第一个已知局面 \(l\) 和右边的第一个已知局面 \(r\) 有关。

记事件 \(X\) 为小 \(P\) 在第 \(x\) 局中获胜，事件 \(L,R\) 为第 \(l,r\) 获胜者和已知相符。则需要求的是：

\[P(X|L,R)=\frac{P(L,R|X)\times P(X)}{P(L,R)}\\
\]

\[=\frac{P(L,R|X)\times P(X)}{P(L,R)}\\
\]

\[=\frac{P(L|X)\times P(R|X)\times P(X)}{P(R|L)\times P(L)}\\
\]

\[=\frac{\frac{P(X|L)\times P(L)}{P(X)}\times P(R|X)\times P(X)}{P(R|L)\times P(L)}\\
\]

\[=\frac{P(X|L)\times P(R|X)}{P(R|L)}\\
\]

我们发现一个区间内的分母是一个定值，因此可以将每个区间内的分母提出，即：

\[\sum_{x=l+1}^{r-1}\frac{P(X|L)\times P(R|X)}{P(R|L)}=\frac{\sum_{x=l+1}^{r-1} P(X|L)P(R|X)} {P(R,L)}
\]

考虑动态更新答案，因此我们每次只需要查询一个区间的如上式子即可。

分母直接对每个点建立矩阵 \(f_i=\left [\begin{matrix}
{1-q_i}&{q_i}\\
{1-p_i}&{p_i}\\
\end{matrix} \right]\) 线段树维护区间矩阵乘积即可。

考虑维护分子，设矩阵 \(g_i=\left [\begin{matrix}
{0}&{q_i}\\
{0}&{p_i}\\
\end{matrix} \right]\) ，发现 \(P(X|L)P(R|X)\) 写成矩阵形式后就是 \((\prod_{i=l+1}^{x-1}\limits f_i)\times g_x \times (\prod_{i=x+1}^{r}\limits f_i)\)

考虑维护 \(\sum_{x=l+1}^{r-1}\limits [(\prod_{i=l+1}^{x-1}\limits f_i)\times g_x \times (\prod_{i=x+1}^{r}\limits f_i)]\)，设 \(F=\prod_{i=l+1}^{r}\limits f_i\)，\(G=\sum_{x=l+1}^{r-1}\limits [(\prod_{i=l+1}^{x-1}\limits f_i)\times g_x \times (\prod_{i=x+1}^{r}\limits f_i)]\)。

对于线段树节点 \(i\)，显然有：

\[F_i=F_{2*i}\times F_{2*i+1}\\
G_i=G_{2*i}\times F_{2*i+1}+F_{2*i}\times G_{2*i+1}
\]

点击查看代码 
#include<bits/stdc++.h>
using namespace std;
int n,m;
char type[3],op[15];
struct mat{
	double a[2][2];
	mat(int _x=1){a[0][0]=a[1][1]=_x;a[0][1]=a[1][0]=0;}
	inline double* operator [](int t){
		return a[t];
	}
	inline void print(){
    	printf("%lf %lf\n",a[0][0],a[0][1]);
    	printf("%lf %lf\n",a[1][0],a[1][1]);
	}
	inline mat operator *(mat b){
		mat res(0);
		for(int i=0;i<2;i++){
			for(int j=0;j<2;j++){
				res[i][j]=(res[i][j]+a[i][0]*b[0][j]);
				res[i][j]=(res[i][j]+a[i][1]*b[1][j]);
			}
		}return res;
	}
	inline mat operator +(mat b){
		mat res(0);
		for(int i=0;i<2;i++){
			for(int j=0;j<2;j++)res[i][j]=a[i][j]+b[i][j];
		}return res;
	}
}f[800005],g[800005];
double p[200005],q[200005];
void build(int l=0,int r=n,int i=1){
	if(l==r){
		f[i][0][0]=1-q[l];f[i][0][1]=q[l];f[i][1][0]=1-p[l];f[i][1][1]=p[l];
		g[i][0][0]=0;g[i][0][1]=q[l];g[i][1][0]=0;g[i][1][1]=p[l];
		return ;
	}
	int mid=(l+r)>>1;
	build(l,mid,i<<1);build(mid+1,r,i<<1|1);
	f[i]=f[i<<1]*f[i<<1|1];
	g[i]=g[i<<1]*f[i<<1|1]+f[i<<1]*g[i<<1|1];
}
mat F,G;
void query(int fr,int to,int l=0,int r=n,int i=1){
	if(fr>r||to<l)return ;
	if(fr<=l&&to>=r){
		G=G*f[i]+F*g[i];F=F*f[i];
		return ;
	}
	int mid=(l+r)>>1;
	query(fr,to,l,mid,i<<1);query(fr,to,mid+1,r,i<<1|1);
}
map<int,bool> mp;
inline double calc(int l,int r,bool op1,bool op2){
	F=mat(1);G=mat(0);
	query(l+1,r);
	return G[op1][op2]/F[op1][op2];
}
double ans;
int main(){
	scanf("%d%d%s",&n,&m,type);
	scanf("%lf",&p[1]);
	for(int i=2;i<=n;i++)scanf("%lf%lf",&p[i],&q[i]);
	p[0]=1;q[0]=1;++n;
	build();mp[0]=1;mp[n]=0;
	ans=calc(0,n,1,0);
	while(m--){
		int i,c;
		scanf("%s",op);
		if(op[0]=='a'){
			scanf("%d%d",&i,&c);
			auto r=mp.upper_bound(i),l=r;l--;
			ans-=calc(l->first,r->first,l->second,r->second);
			ans+=calc(l->first,i,l->second,c);ans+=calc(i,r->first,c,r->second);mp[i]=c;
		}
		else {
			scanf("%d",&i);
			auto it=mp.find(i),l=it,r=it;l--;r++;
			ans-=calc(l->first,it->first,l->second,it->second);ans-=calc(it->first,r->first,it->second,r->second);
			ans+=calc(l->first,r->first,l->second,r->second);mp.erase(it);
		}printf("%lf\n",ans);
	}

	return 0;
}

方差相关

方差等于平方的期望减去期望的平方

\[Var(x)=\frac{1}{n}\sum_{i=1}^{n}(x_i-\overline{x})^2\\
=\frac{1}{n}\sum_{i=1}^{n}(x_i^2-2\overline{x}*x_i+\overline{x}^2)
\]

\[=\frac{\sum_{i=1}^{n} x_i^2}{n}-2\overline{x}*\frac{\sum_{i=1}^{n} x_i}{n}+\frac{\sum_{i=1}^{n} \overline{x}^2}{n}
\]

\[=\frac{\sum_{i=1}^{n} x_i^2}{n}-2\overline{x}^2+\overline{x}^2
\]

\[=E(x^2)-E(x)^2
\]

矩形覆盖

给定一个大小为 \(n*m\) 的矩形，某一些格子上有物品，共有 \(k\) 个物品，现在等概率选一个子矩形，求子矩形内物品个数的方差。

考虑 \(Var(x)=E(x^2)-E(x)^2\)，分别求 \(E(x)\) 和 \(E(x^2)\)。

\(E(x)\) 较为好求，即统计每一个点的贡献即可。

求 \(E(x^2)\) 可以考虑对于每一对点，统计包含它们的矩形的个数，考虑分类讨论一下它们是左上右下的关系还是右上左下的关系，扫描线计算一下即可。

【UER #6】逃跑

还是考虑求 \(E(x^2)\) 和 \(E(x)\) ：

预处理 \(g[t][x][y]\) 表示 \(t\) 时刻从 \((0,0)\) 到达点 \((x,y)\) 的方案数。

设 \(f[i]\) 为 \(i\) 时刻期望经过的的点数，有：

\[f[i]=(w1+w2+w3+w4)^i-\sum_{j=0}^{i-1}f[j]\times g[i-j][0][0]
\]

则 ：

\[E(x)=\sum_{i=0}^{n}f[i]\times (w1+w2+w3+w4)^{n-i}
\]

考虑统计同时经过的点对数的方式求出 \(E(x^2)\) ：

设 \(h[t][x][y]\) 表示在 \(t\) 时间内经过了 \((a,b)\) 并在 \(t\) 时刻到达 \((a+x,b+y)\) 的方案数，有：

\[h[t][x][y]=\sum_{i=0}^{t-1}f(i)\times g[t-i][x][y]-h[i][-x][-y]\times g[t-i][x][y]-h[i][x][y]\times g[t-i][0][0]
\]

由于我们统计的是点对，因此 \(\sum_{t=0}^{n}\limits \sum_{x=-n}^{n}\limits \sum_{y=-n}^{n}\limits h[t][x][y]\times (w1+w2+w3+w4)^{n-t}\) 求出的实际上是 \(E(\frac{n\times (n-1)}{2})\) 。

因此 \(E(x^2)=E(x)+2\times \sum_{t=0}^{n}\limits \sum_{x=-n}^{n}\limits \sum_{y=-n}^{n}\limits h[t][x][y]\times (w1+w2+w3+w4)^{n-t}\) 。

题目要求 \(Var(x)\times (w1+w2+w3+w4)^n\) ，而 \(\text{std}\) 求的是 \(Var(x)\times (w1+w2+w3+w4)^{2\times n}\)。

因此最终输出应该输出 \(E(x^2)\times (w1+w2+w3+w4)^n-E(x)^2\) 。

点击查看代码 
#include<bits/stdc++.h>
using namespace std;
int n;
long long w1,w2,w3,w4,sum;
const long long md=998244353;
int g[105][205][205],h[105][205][205],f[105],pwr[105];
long long E1,E2;
int main(){
	scanf("%d",&n);
	scanf("%lld%lld%lld%lld",&w1,&w2,&w3,&w4);sum=(w1+w2+w3+w4)%md;
	g[0][n][n]=1;
	for(int t=0;t<n;t++){
		for(int i=0;i<=2*n;i++){
			for(int j=0;j<=2*n;j++){
				if(!g[t][i][j])continue;
	//			cout<<t<<" "<<i<<" "<<j<<" "<<g[t][i][j]<<endl;
				g[t+1][i-1][j]=(g[t+1][i-1][j]+w1*g[t][i][j])%md;
				g[t+1][i+1][j]=(g[t+1][i+1][j]+w2*g[t][i][j])%md;
				g[t+1][i][j-1]=(g[t+1][i][j-1]+w3*g[t][i][j])%md;
				g[t+1][i][j+1]=(g[t+1][i][j+1]+w4*g[t][i][j])%md;
			}
		}
	}
	pwr[0]=1;
	for(int i=1;i<=n;i++)pwr[i]=pwr[i-1]*sum%md;
	for(int i=0;i<=n;i++){
		f[i]=pwr[i];
		for(int j=1;j<=i;j++)f[i]=(f[i]-1ll*f[i-j]*g[j][n][n])%md;
	//	cout<<i<<" "<<f[i]<<endl;
		E1=(E1+1ll*f[i]*pwr[n-i])%md;
	}E2=E1;
	for(int t=0;t<=n;t++){
		for(int i=0;i<=2*n;i++){
			for(int j=0;j<=2*n;j++){
				if(i==n&&j==n)continue;
				for(int k=0;k<t;k++)h[t][i][j]=(h[t][i][j]+1ll*f[k]*g[t-k][i][j]-1ll*h[k][2*n-i][2*n-j]*g[t-k][i][j]-1ll*h[k][i][j]*g[t-k][n][n])%md;
				E2=(E2+2ll*h[t][i][j]*pwr[n-t])%md;
			}
		}
	}
	printf("%lld",((E2*pwr[n]-E1*E1)%md+md)%md);

	return 0;
}

概率生成函数

对于任意取值在非负整数集上的离散随机变量 \(x\)，它的概率生成函数为：

\[F(z)=\sum_{i=0}^{\infty}P(x=i)z^i
\]

概率生成函数的一些性质

1.

\[F(1)=\sum_{i=0}^{\infty}P(x=i)=1
\]

2.

\[F'(1)=\sum_{i=0}^{\infty}P(x=i)*i=E(x)
\]

3.

\[F^{(k)}(1)=\sum_{i=0}^{\infty}P(x=i)*i^{\underline{k}}=E(x^{\underline{k}})
\]

4.

\[E(x^{k})=E(\sum_{i=0}^{k}\left \lbrace\frac{k}{i}\right\rbrace x^{\underline{i}})=\sum_{i=0}^{k}\left \lbrace\frac{k}{i}\right\rbrace E(x^{\underline{i}})=\sum_{i=0}^{k}\left \lbrace\frac{k}{i}\right\rbrace F^{(i)}(1)
\]

[CTSC2006]歌唱王国

设 \(f[i]\) 表示第 \(i\) 轮结束的概率，\(g[i]\) 表示到了第 \(i\) 轮还没有的概率。

设生成函数：

\[F(x)=\sum_{i=0}^{\infty}\limits f[i]\cdot x^i，G(x)=\sum_{i=0}^{\infty}\limits g[i]\cdot x^i
\]

显然：

\[f[i]=g[i-1]-g[i]
\]

有：

\[F(x)=x\cdot G(x)-G(x)
\]

我们考虑在一个长度为 \(i\) 且没有结束的序列后面生成一个目标序列，直接在序列后面生成的概率是 \(g[i]\times (\frac{1}{n})^m\) 。

假如目标序列的 \(a_{1\to j}=a_{m-j+1\to m}\)，那么我们在长度为 \(i+j\) 且结束了的序列后面插入 \(m-j\) 个字符也可以生成相同的序列，概率为 \(\sum_{j=1}^{m}\limits [a_{1\to j}=a_{m-j+1\to m}]f[i+j]\times (\frac{1}{n})^{m-j}\)

记 \(kmp[i]=[a_{1\to i}=a_{m-i+1\to m}]\) 。

有：

\[g[i]\times (\frac{1}{n})^m=\sum_{j=1}^{m}\limits kmp[j]\times f[i+j]\times (\frac{1}{n})^{m-j}
\]

进而：

\[G(x)\cdot \left(\frac{x}{n}\right)^m=\sum_{j=1}^{m} kmp[j]\times F(x)\cdot \left(\frac{x}{n}\right)^{m-j}
\]

对 \(F(x)=x\cdot G(x)+G(x)\) 求导可得：

\[F'(x)=G(x)+x\cdot G'(x)-G'(x)
\]

令 \(x=1\)，有：

\[F'(1)=G(1)+G'(1)-G'(1)=G(1)
\]

联合第二个方程，有：

\[F'(1)\cdot \left(\frac{1}{n}\right)^m=G(1)\cdot \left(\frac{1}{n}\right)^m=\sum_{j=1}^{m} kmp[j]\cdot F(1)\cdot \left(\frac{1}{n}\right)^{m-j}\\
F'(1)=\sum_{j=1}^{m} kmp[j]\cdot n^j\\
\]

点击查看代码 
#include<bits/stdc++.h>
using namespace std;
int n,t,m;
int a[100005],nxt[100005],pwr[100005];
const int md=10000;
int main(){
	scanf("%d%d",&n,&t);
	pwr[0]=1;
	for(int i=1;i<=1e5;i++)pwr[i]=pwr[i-1]*n%md;
	while(t--){
		scanf("%d",&m);
		for(int i=1;i<=m;i++)scanf("%d",&a[i]);
		for(int i=2,j=0;i<=m;i++){
			while(j&&a[i]!=a[j+1])j=nxt[j];
			if(a[i]==a[j+1])j++;
			nxt[i]=j;
		}int ans=0;
		for(int i=m;i;i=nxt[i])ans=(ans+pwr[i])%md;
		if(ans<1000)printf("0");
		if(ans<100)printf("0");
		if(ans<10)printf("0");
		printf("%d\n",ans);
	}

	return 0;
}

Dice

有一个 \(m\) 面的骰子，求扔连续 \(n\) 次相同就结束的期望步数和扔连续 \(n\) 次结果不同就结束的期望步数。

对于第一问，可以列出：

\[\begin{cases}
F(x)=x\cdot G(x)-G(x)\\
\\
G(x)\cdot (\frac{x}{m})^{n-1}=\sum_{i=1}^{n}\limits F(x)\cdot (\frac{x}{m})^{n-i}
\end{cases}
\]

有：

\[F'(1)=G(1)\to F'(1)\cdot \left(\frac{1}{m}\right)^{n-1}=\sum_{i=1}^{n}\limits F(1)\cdot \left(\frac{1}{m}\right)^{n-i}
\]

解得：

\[F'(1)=\sum_{i=1}^{n}\limits m^{i-1}=\frac{m^n-1}{m-1}
\]

第二问：

\[\begin{cases}
F(x)=x\cdot G(x)-G(x)\\
\\
G(x)\cdot (\frac{x}{m})^{n}\cdot m^{\underline{n}}=\sum_{i=1}^{n}\limits F(x)\cdot (\frac{x}{m})^{n-i}\cdot (m-i)^{\underline{n-i}}
\end{cases}
\]

有：

\[F'(1)=G(1)\to F'(1)\cdot(\frac{1}{m})^{n}\cdot m^{\underline{n}}=\sum_{i=1}^{n}\limits F(1)\cdot (\frac{1}{m})^{n-i}\cdot (m-i)^{\underline{n-i}}
\]

解得：

\[F'(1)=\sum_{i=1}^{n}\limits \frac{m^{i}}{m^{\underline{i}}}
\]

例题

CF1187F Expected Square Beauty

点击查看代码 
#include<bits/stdc++.h>
using namespace std;
int n;
int l[200005],r[200005];
const long long md=1e9+7;
long long R[200005],pre[200005],suf[200005];
inline long long pwr(long long x,long long y){
	long long res=1;
	while(y){
		if(y&1)res=res*x%md;
		x=x*x%md;y>>=1;
	}return res;
}
long long ans;
int main(){
	scanf("%d",&n);
	for(int i=1;i<=n;i++)scanf("%d",&l[i]);
	for(int i=1;i<=n;i++)scanf("%d",&r[i]);
	for(int i=0;i<n;i++){
		int len=max(0,min(r[i],r[i+1])-max(l[i],l[i+1])+1);
		R[i]=(1-len*pwr(r[i]-l[i]+1,md-2)%md*pwr(r[i+1]-l[i+1]+1,md-2))%md;
	}
	for(int i=0;i<n;i++)pre[i]=(pre[i-1]+R[i])%md;
	for(int i=n-1;~i;i--)suf[i]=(suf[i+1]+R[i])%md;
	for(int i=2;i<n;i++)ans=(ans+R[i]*pre[i-2])%md;
	for(int i=0;i<n;i++)ans=(ans+R[i]*suf[i+2])%md;
	for(int i=0;i<n;i++)ans=(ans+R[i])%md;
	for(int i=0;i+1<n;i++){
		int len=max(0,min(r[i],min(r[i+1],r[i+2]))-max(l[i],max(l[i+1],l[i+2]))+1);
		ans=(ans+2*(R[i]+R[i+1]-1+len*pwr(r[i]-l[i]+1,md-2)%md*pwr(r[i+1]-l[i+1]+1,md-2)%md*pwr(r[i+2]-l[i+2]+1,md-2)%md)%md)%md;
	}
	printf("%lld",(ans+md)%md);

	return 0;
}

原文链接: https://www.cnblogs.com/A-Quark/p/15834749.html

欢迎关注

微信关注下方公众号，第一时间获取干货硬货；公众号内回复【pdf】免费获取数百本计算机经典书籍

    概率与期望

原创文章受到原创版权保护。转载请注明出处：https://www.ccppcoding.com/archives/187488

非原创文章文中已经注明原地址，如有侵权，联系删除

关注公众号【高性能架构探索】，第一时间获取最新文章

转载文章受原作者版权保护。转载请注明原作者出处！

赞 (0)
adminadmin
0 0
上一篇 2023年2月12日 上午11:18
下一篇 2023年2月12日 上午11:18

相关推荐