第二章 蒙特卡洛方法
.0 概率与统计 古典概率 : 在相同的实验条件下, 随机事件 A,B 按各自确定的概率发生 - 和概率 A.OR.B : P(A+B) = P(A) + P(B) P(AB) - 与概率 A.AND.B: P(A*B) = P(A B)*P(B) = P(B A)*P(A) - 条件概率 P(A B) = 在随机事件 B 发生的条件下,A 发生的概率 - 互斥 P(A*B) = 0,ie 随机事件 AB 不能在同一实验中同时发生 - 相互独立 P(A*B) = P(A)*P(B),ie P(A)=P(A B)=P(A )
全概率公式 : 随机事件 A 构成互斥完备集合 {A i }, 则任意随机事件 B 可表述为 P(B)= P(B A i)p(a i) i 贝叶斯 Bayes 公式 : P(A B ) = i P(AiB) P(B A i)p(a i) = P(B) P(B A j)p(a j) j 3
随机变量 X: ) 在相同的确定实验条件下, 对 X 的观测无法给出单一固定值 ; ) 必须依据遍举测量原则, 对所有可能取值给出发生概率 离散变量举例 : 3MeV 光子入射屏蔽铅板的全吸收反应过程 反应类型 X : 光电效应 Compton 散射电子对产生 反应几率 : ε ε ε3 其中 i.e. ε + ε + ε3 = 00% n { }, = i i = i = X p 4
连续型随机变量 : X 在连续区间取值, 其取某确定值的概率由分布密度函数给出 f ( d ) = P [ < X< + d] 分布函数 F ( ) = f ( ) d 则有 + F( ) f d + = ( ) 5
联合分布密度 : 相互独立 : 描述两个 (i.e. 多维 ) 随机变量 X 与 Y 的相互关联 ( ) + f = f (, y) dy ( ) = f, y f ( ) f ( y) 6
函数的分布密度 : 随机变量 X 密度函数 f(), 其函数 Y=Y(X) 的密度函数 f ( ) d = g( y) dy 几率密度相同 则 ( ) f( ( y)) g y = y 变量变换 Jaccobi 7
随机变量的特征值 ) 期望值 (mean): 出现几率最大或概率中心的观测值 + Ey [ ( )] yˆ = y( ) f ( d ) ) 方差 (standard deviation): 随机变量 分布对期望值的离散程度 + Dy [ ( )] σ ( [ ( )]) ( ) [( ˆ y = y Ey f d= E y y) ] 3) 特征运算 : E c c D c E c ce D c c D [] =, [] = 0, [ ] = [], [ ] = [] E [ + ] = E [ ] + E [ ] D [ + ] = D [ ] + D [ ] + Cov [, ] 8
几种著名分布 ) 二项式分布 (Binominal): 发生几率为 p, 不发生为 q=(-p), 则 N 次试验中出现 k 次的几率 N! k Pk [ ] = p p k! N k! 其中 k=0,,,3, 0 p, p+q ( ) ( ) N k E[ k] = Np, D[ k] = Np( p) 例 : 反应触发率 (trigger rate) 定义为 ε=k/n, 求其期望值 E[ε] 与方差 D[ε] 9
) 泊松分布 (Poisson): 在相同实验条件下, 相同时间内, 随机过程发生 k 次的几率 Pk [ ] 其中关于分布参数 λ 有 λ k λ = e k 0, λ > 0 k! E[ k] = D[ k] = λ 当 λ 时,Poisson 分布过渡到 Gaussian 分布 0
3) 高斯分布 (Gaussian): ( μ) σ f() = e N( X; μ, σ) < <+ πσ 标准化 (-μ)/σ, 则正态分布 f() = e N(0,) π E ( ) = μ D ( ) = σ
分布函数 对称分布 y ( y) e d π Φ = Φ ( y) = Φ( y) 则, 如 [a,b] 对 μ 对称, 有 b μ a μ b μ Pa ( b) =Φ( ) Φ ( ) = Φ( ) σ σ σ
Gaussian 计数 Pk [ ] = λ k λ e k! E(k)=D(k)=λ λ e.g. λ>0 N(k;m,m) μ=m=λ σ =m=λ 则, 当统计计数时,N (N 0), 过渡至高斯分布 计数期望值 μ=ν 均方根方差 σ= N 3
4) 指数分布 (Eponential): 描述自由粒子寿命, 或粒子平均自由程分布函数 分布函数 f ( ) e λ / λ = > 0, λ > 0 z F ( z) = f ( ) d = e z / λ 5) 均匀分布 (uniform): 其中 f ( ) =, a b b a a + b ( b a) E ( ) =, D ( ) = 4
大数法则 : 在函数 f() 定义域 [a,b] 内, 以均匀概率分布随机地取 N 个数 i, 函数值之和的算术平均收敛于函数的期望值 lim N N b f ( ) lim ( ) i IN f d I N N b a a i= 在抽取足够多的随机样本后, 积分的蒙特卡洛估计值 ( 左边 ) 将收敛于该积分的正确结果 ( 右边 ) 即随机变量统计量为 b N b a E[ y ( )] = y ( ) f( d ) lim y ( i) f( i) a N N i= 5
中心极限定理 : 大量微弱因素累加而成的随机变量服从单一正态分布 例 :n 个相互独立分布各异的随机变量,n, 则总和服从正态分布 N μ, σ ep μ /σ σ π ( ) = ( ) Gaussian 分布随机测量报道 置信水平 = μ ± k σ μ μ P{ μ kσ μ+ kσ} = P k k = P k σ σ = Φ( k) Φ( k) = Φ( k) 6
统计量 : 随机变量 X(μ,σ) 一组测量样本 { i } 的函数 y = y({ i }) 例, 样本平均值 N i y = N i = 作为 N 维相互独立 ( 测量 ) 随机变量的函数,y 亦为随机变量, 亦存在分布 E ( ) = E ( ) = μ D ( ) σ D ( ) = = N N 收敛于期望值 ( 大数定理 ) 期望值测量误差 ( 中心极限 ) 7
. Monte Carlo 方法 针对待求问题, 根据物理现象本身的统计规律, 或人为构造一合适的依赖随机变量的概率模型, 使某些随机变量的统计量为待求问题的解, 进行大统计量 N 的统计实验方法或计算机随机模拟方法 理论依据 : - 均匀分布的算术平均收敛于真值 ( 大数法则 ) - 置信水平下的统计误差 ( 中心极限 ) 待求问题 : ) 自然界中本身存在的随机过程, 如粒子衰变过程 粒子在介质中的输运过程等 ) 以慨率模型来解决不直接具有随机性的确定性问题, 如求 π 求积分 8
例. Buffon 投针实验求 π(777 年 ) ) 平行线间距 = 针长 =s;) 针与平行线垂线方向夹角 a 则相交概率为 s cos a P = = s cos a 3) 各项同性均匀抛针,i.e. 夹角 a 在 [0,π] 均匀分布 π E( cos a) = cos a f( a) da= 0 π 4) 设 N 次抛针,M 次相交, 则相交概率的期望值为 E(cos a) = π M N π = N M (N 大数定理 ) 9
问 :π 的测量精度? 服从二项式分布, 单次相交概率 则 C p ( p) M M N M N EM ( ) = Np M p = = N π DM Np p σ M ( ) = ( ) ( ) ( ) ( ) ( M ) ( ) p σ p = σ = σ M = p N N N 0
Qπ = p d 4 ( ) ( ) ( ) ( π p p σ π ) σ ( p) π π = = = 4 σ π = dp p N N.37 N N π = ±.37 M N 对真值的测量精度与测量次数平方根反比, 即 0 6 次实验才精确到 0-3
例. 投点法求定积分 I = 0 f ( ) d 随机地向 [0,], y [0,y ma ] 正方形内投点 N, 统计落在曲线下的点数 M, 当总掷点数 N 时 N-M M I = M N 建立恰当的概率模型, 即确定某个随机事件 A 或随机变量 X, 使得待求问题的解等于随机事件出现的概率或随机变量的数学期望值 然后重复进行多次的随机实验, 对结果进行统计平均, 求出 A 出现的频数或 X 的平均值作为问题的近似解 这种方法也叫做蒙特卡洛模拟
. 伪随机数 进行计算机模拟需要大样本的均匀分布随机数数列, 如何获得? - 真随机数 : 由随机物理过程来产生, 例如 : 放射性衰变 电子设备的热噪音 宇宙射线的触发时间等等 - 伪随机数 : 由计算机按递推公式大量产生 冯. 诺曼平方取中法 乘加同余 n+ = ( a ξ n r ( ) r n + n mod n = = ξ = / n / m n + c)(mod n r m) ξ [0,] n 3
伪随机数的统计检验 均匀性检验 : - [0,] 分成 k 个相等子区间, 进行 N 次抽样, 投入各子区间 - 如均匀, 则各区间落入数 N i 应为 N = N i k i= N = = = μ i σ N k 计数 Poisson Gaussian -N i 可视为 (μ,σ) 的一组无关样本测量, 服从 χ ( N μ) k i = χ i= σ ( k) 则 χ k 4
独立性检验 : 即 ξ i 与 ξ i+ 的前后无关性 - [0,] 上进行 N 次抽样, 分成两个序列 X : ξ, ξ,... ξ,..., ξ 3 i N Y : ξ, ξ,... ξ,..., ξ 4 i N - 在 XY 平面内划分 k k 方格, 如独立, 则各格内落入数应为 - 则服从 χ 分布 χ N = n ij k i, j= n ij = = = μ σ N k ( n μ) k ij = χ k ij= σ ( ) 满足以上统计性检验的递推抽样序列, 可视为 [0,] 均匀分布伪随机数 5
.3 任意分布的伪随机变量的抽样 A. 离散型分布 随机变量 本证值 本证值 概率密度 X :{,,,..., } 3 f :{ p, p, p,..., p } 3 N N 归一化, 几率守恒 分布函数 N i= p i = F( ) = pi 0 i 则取 [0,] 均匀分布随机数 ξ, 按分布函数不等式抽样 η ( ) ξ F( ) F < η = j j j 6
离散变量举例 : 3MeV 光子入射屏蔽铅板的全吸收反应过程 反应类型 X : 光电效应 Compton 散射 电子对产生 反应截面 : σ σ σ3 反应几率 : ε= σ/σ ε= σ/σ ε3= σ3/σ 其中,σ+σ+σ3=σ,ε + ε + ε3 = 00% ξ ξ < ε? N ξ < ε? N 电子对产生 Y Y 光电效应 Compton 7
B. 连续型分布. 反函数直接抽样法 : 类比于离散型随机抽样, 求分布函数 0 F ( ) = f ( ) d 则令 ξ = F( η) η = F ( ξ ) 例对指数分布的直接抽样 f ( ) = λ λ e, > 0, λ > 0, 其它. 0 8
例对指数分布的直接抽样 积分得到分布函数 F f ( ) λ λ e, > 0, λ > 0 = 0, 其它. λt λ ( ) f ( t) dt = λe dt = e = 0 令 λη ( η) = ξ = F e 则指数分布的随机变量抽样为 η = ln ( ξ) lnξ λ λ 9
. 变换抽样法 : 随机变量 按 f() 抽样已知, 其函数 y=y() 的概率密度有 f ( d ) = g( ydy ) g( y) = f( ( y)) d dy 二维情况 : 联合分布密度函数 g(u,v) 的随机变量 (η,δ) 抽样方法已知, 则其函数 = g ( u v) y = g ( u, v), 抽样可依据 f(, y) = g( u(, y), v(, y)) u v u y v y 30
例 Gaussian 正态分布 f f ( ) 可由标准正态分布代替 ep π σ ( μ) = σ ( ) = ep{ / } π η = σδ + μ 解 : ) 相互独立均匀抽样 u,v, i.e. g(u,v)== ) 构造二维随机变量 ( π ) ( π ) = lnucos v y = lnusin v ( u ep y ) = + v = tan / π ( y ) 3
则,y 联合分布密度函数为 ( ) f, y = g u, y, v, y J ( ) ( ) ( ) 代入 Jacobi 行列式, 有 其中, f y y f f y π ( ) (, = ep + ) = ( ) ( ) f ( ) = ep{ / } π f ( y) = ep{ y / } π 3
改进的 Maraglia 方法 : () 产生 [0,] 区间上的独立均匀分布随机数 u,v u,v [0,] () 计算 w ( u ) + ( ) = v ω=(u-) +(v-) (3) 如果 ω>, 回到步骤 (); 否则, 执行 (4) (4) 计算 z = [ ln( w) / w] / (5) 得到二维相互独立标准正态分布抽样 = uz, y = vz ω>? Y z=[-lnω/ω] / X=uz N 舍选效率 ~% 33
3. 舍选抽样法 : 第三类舍选 不易直接抽取随机变量 X 分布 f(), 而与随机变量 Y 联合分布密度 g(,y) 存在 + h( ) f ( ) g(, y) dy L g(, y) dy (.3.38) 联合分布 Page5 数学上可构造 其中 函数 h() 取值在 y 定义域上 常数 L 的定义为 L + h( ) = g(, y) ddy > 34
则有 { η ( )} ηy η p h = h( t ) { η } { η ( )}, ηy η p{ η h( η )} p h dt g( t, t) dt h( t ) = = L g( t, t ) dt dt + h( t ) dt g( t, t ) dt p Page Bayes 公式 y 子样 η 概率等于在 η y h(η ) 条件下,η 出现的概率 35
第三类舍选抽样步骤 : ⑴ 由联合分布密度函数 g(,y) 抽取随机向量 η, η ) ( y 按 g(,y) 抽取 (η,η y ) ( ) ⑵ 判别 η y h η 是否成立 ; 若不成立, 返回 ⑴ ⑶ 取分布密度函数 f() 的抽样值 η = η η y h(η )? Y N 抽样效率 /L η=η 36
例 : 各向同性方位角余弦的抽样 解 : 物理含义 φ [0,π] 各向同性 η=cos(πξ) 密度函数 f( ) = π 0 < 其它 cos F( ) = f ( t) dt = dt π = t π 需要计算 cos 函数 直接抽样 cos η ξ = η = cos( πξ ) π 37
解 3 独立抽取均匀分布 ξ 与 ξ, 定义 y = ξ ξ = ξ + ξ + ξ ξ ξ = ( ) y + ξ = ( ) y (ξ,ξ ) 到 (,y) 变换抽样 gy (, ) = ϖξ (, ξ) J = J = 则构造 (,y) 联合密度函数 gy (, ) = 4 4, <,0< y< 0, 其它 38
类比构造 h() 与 L h( ) f( ) = L g(, y) dy π h( ) h( ) = L dy = L dy 4 4 0 0 则有 L 4 =, h( ) = π 39
抽样 : 改进 : ξ, ξ ξ, ξ y = = ξ ξ ξ + ξ + ξ ξ N y = ξ = ξ N y h( ) = η = η = ξ Y 效率 /L=π/4 ξ ξ + ξ η = + y Y y + y y η = + y 40
第二类舍选 判据函数 h( ) f ( ) L g(, y) dy 抽样效率 如,y 相互独立, ( y) g ( ) g ( y) g, =, 且 y 为均匀分布 抽样密度函数 则有 g ( y), y [0,] = 0, 其它 h ( ) ( ) ( ) 在 y 的定义域上 ( ) 0 h f = Lg g ( y) dy = L h( ) g ( ) 效率判据抽样 4
则对任意密度函数构造 ( ) ( ) f f ( ) = L g( ) L h( ) g( ), Lg 效率判据抽样 h ( ) f ( ) ( ) Lg Y 均匀分布 g (y)= 0 h(), 即 L = f ( ) ma h ( ) > 定义域 ξ, 按 g() 抽取 η g 密度概率函数 g(), 易于抽样 y=ξ h(η g )? N Y 抽样效率 /L η =η g 4
例采用第二类舍选抽样标准正态分布 f ( ) = ep π ( < < + ) 解 : 不存在反函数, 凑二类舍选函数形式 指数直接抽样 : g( ) e ηg = lnξ 判据 : 效率 : 偶函数 : L e ; π ( ) h ( ) ep{ } ( 0 < < + ) ( ηg ) ln η Y = η g ξ N 43
第一类舍选 二类舍选中, 如抽样概率密度函数 g() 均匀分布, 则 L = f ( ) = ma 定义域 [ a, b] λ ξ,ξ δ = a+ ( b a) ξ λf ( ) ξ λf Y ( ) δ N η = δ 44
例对随机变量抽样 f ( ) =,0, 0, 其它. ( ) 解 : 反函数法 F =, 直接抽样 = ξ 第一类舍选 ξ,ξ 3 改进 η = ma( ξ, ξ ) λ f ( δ ) = ξ ξ ξ η Y = ξ N 4 推广 f ( ) n n, [0,], n =,,... = 0, 其它. ( ξ, ξ,..., ) = ma ξ n 45
4. 复合抽样法 : 复合分布密度函数 + f ( ) = g ( y) h( y) dy 根据密度函数 h(y) 抽样 y h 根据条件密度函数 g( yh) 抽样 f g = ξ g ( y h ) 46
加分布抽样 复概率密度函数 f ( ) = hn( ) p n n 其中 抽样步骤 : 0< p n <, p n = ; hn ( ) = n n i= + 取 [0,] 区间上均匀分布随机数 ξ n p i < ξ p i=, 解下不等式求 n i 找到对应的 h n (), 对其抽样得到 η = η h n 47
例 : 球壳均匀分布的抽样 f ( 3 r r ) =, R R R 3 0 r R 3 0 例 : 积分求分布函数, 直接抽样 [( ) ] 3 3 3 R R + / 3 η = ξ 0 R0 变量代换 r = ( R R ) + R, 0 0 λ = R + + RR0 R0 密度函数 f ( ) = ( R λ R 0 ) 3 + 3R 0 ( R λ R 0 ) + 3R λ 0 48
49
减分布抽样 复概率密度函数 f( ) = Ag ( ) A g ( ), A, A > 0 取 m = g ( ) min, [ a, b] g ( ) g (),g ()>0 概率密度函数 则 即 g ( ) 0 < f ( ) = g( ) A A g( )( A Am) g( ) f( ) A A g ( ) 0 < h ( ) = ( A Am) g( ) A Am A Am g( ) 50
亦即第二类舍选 f ( ) = ( A A m) h ( ) g ( ) 效率 L >? (0,], 判据概率密度函数, 抽样 另第二类舍选 其中 A A m f h g m ( ) = ( ) ( ) < Am g ( ) Am h = 0 ( ) A Am g( ) A Am 5
乘加分布抽样 考虑 n= 情况 f ( ) = H n( ) gn( ) n f ( ) = H ( ) g ( ) + H ( ) g ( ) (.3.69) 令 b p = H ( ) g ( ) d, p = H( ) g( ) d, p + p a b a = 则加分布 H( ) H( ) f ( ) = p g( ) + p g( ) = pg ( ) + pg ( ) p p 问题 : 需要计算 p 与 p 积分 5
Give a new try 其中 M H( ) M H ( ) f ( ) = ( M + M ) g( ) + g( M + M M M + M M M i = ma [ a, b] H () i ) 令 取 p M =, p = M M + M M+ M L + = L = M M Hi( ) 0 < hi ( ) = M i 则有 [ ] [ ] f ( ) = p L h( ) g ( ) + p L h ( ) g ( ) 加分布 + 第二类舍选 53
乘减分布抽样 f ( ) = H ( ) g ( ) H ( ) g ( ), [ a, b] 令 m H ( ) g( ) min [, b] H ( ) g ( ) = M = a [ a, b] ma H ( ) 有 H ( ) g( ) 0 < f ( ) = H( ) g( ) H( ) g( )( m) M( m) g( ) H( ) g( ) 再令 H ( ) g ( ) 0 < h( ) = H( ) M( m) g( ) 给出 f ( ) = M ( m) h( ) g ( ) 第二类舍选 54
4. 特殊抽样方法 : 反函数近似 - 问题 : 根据分布函数 F () 直接抽样最为便捷, 但其反函数 F ( y) 难以解析求出 - 解答 : 分布函数 () 性质 近似替代 F F ( y) Q( y) y [0,] lim F ( y 0 y) = lim F ( y y) + 最小二乘拟合 F ( y) Q( y) = a + by + cy + α ( y) ln y + βy ln( y) (.3.86) y 0, y, 55
例 : 标准正态分布 N(0,) y = F ( y) Q( y) e π d e π d 将 [0,] 区间 00 等份 ( k =,,...,99) y k = Q( y) e π d, Q( y k ) y e π d 逐步回归法计算 a = 0.868, b =.6736, c = 0, α = 0.335, β = 0.335 56
极限近似法期望快速正态分布随机变量抽样 ; 利用均匀分布伪随机数 ξ [0,] E ( ξ) =, D( ξ) = (Page4) 则定义 Rn = ξ + ξ +... + ξ n, 则 + + + n n E ( Rn) = n f ( ) d = n d =, D ( Rn) = 据中心极限定理 (Page6-7), 有 n n n DR ( n) n NR ( n;, ): ER ( n) = ER ( n) =, DR ( n) = = n 57
引入服从标准正态分布随机变量 δ n = R n n N( δn;0,) n 通常取 n=, 则 δ = R 6 注意 : 以上 n= 处理仅能描述 < 6 定义域范围内的正态分布抽样 58
5. 多维随机向量的抽样方法 : 舍选法 多维相互关联随机变量的联合分布密度函数 f (,,, n { a b, a b,, a b } n n ), 定义域为 n 如存在有限上限 L = sup f (,,, n ) < + 则可采取第一类舍选 ξ < + f b a ξ a b a ξ a b a ξ a L + + + (( ),( ),,( ) ) n n n n n 抽样效率 E = n L i= ( b i a i ) 59
条件密度法 三维随机向量的联合分布密度函数可写为 f (,, 3) = f( ) f ( ) f3( 3, ) 其中 Bayes 公式, 有 + + f( ) = f(,, 3) dd3 + f( ) = f(,, 3) d3/ f( ) f3( 3, ) = f(,, 3)/[ f( ) f( )] 60
抽样步骤 f ( ) () 由为分布密度函数产生抽样值 η = η = η = () 在 条件下, 由分布密度函数 f ( η ) 产生抽样值 η=, η= f3( 3 η, η) η = (3) 在条件下, 由分布密度函数产生抽样值 3 3 可进一步推广至 n 维随机向量 f (,,, n ) = f( ) f ( ) f3( 3, ) f n ( n,, n ) 得到随机变量 (,,,, ) 3 n T 的一组抽样向量 r η = η η η η n (,,,, ) 3 T 6
例中子入射角服从联合分布密度函数 f ( ϕ, θ ) = ( + 3sinϕ sinθ )sinϕ sin θ, α π / > ϕ 0,0 θ π /, α = (3 + 3) π /. 求余弦抽样 = cos ϕ, y= cosθ 解 : 变量代换 (ϕ,θ) (,y), 联合密度函数为 f y y y f f y (3+ 3) π (, ) = + 3 ( )( ) = ( ) ( ) 则有 + f 3 ( ) f (, y ) dy π (3 3) π 4 3 = = + + ( ) + 3 ( )( ) f ( y ) = y y π 3 + 4 3 6
先对 f ( ) 抽样, 加分布 =η 产生随机数 ξ f( ) = p + 3 p =, 3 + 3 p 4 π p = 3 3 + 3 ξ < 3/(3+ 3)? No 产生随机数 ξ, ξ3 Yes ξ + ξ3 <? No Yes η = ξ 63
再按 f ( y = η) 加分布抽样 y=δ π 4 f 4 ( y = η) = y π 3 π + η 4 3 3 η 3 3 ( ) + y π 3 + η 4 3 ξ < 4 产生随机数 ξ4 π /4 π 3 + ξ 4 3 产生随机数 ξ5, ξ6? No 产生随机数 ξ5, ξ6 ξ 5 + ξ6 <? No No ξ 5 + ξ6 <? Yes Yes δ = ξ 5 64
高维正态分布随机向量抽样 相互独立 f( y, y,, yn) = ep y / ep y / ep y / π π π 可进行 n 次一维变量正态分布抽样 r { } { } { } n η yi, 获得抽样样本向量 (,,,, ) η = η η η η y y y y y 3 n T 65
相互关联 f r (,, n ) = M r n T ( π ) ( χ μ) ( χ μ) M ep r r 其中 r χ = ( ),,, n, T r r μ = E η = μ μ μ { } (,,, ) n T 正定对称 n 阶协方差矩阵 M 的矩阵元 {( ηi μi )( η j μ j )} ji σ ij = E = σ 存在非奇异的下三角矩阵 A M = T AA 则关联抽样可由无关抽样表示 r r r η = μ + Aη y 66
例二维正态分布的抽样 解设 r η = ( η, η ) 其协方差矩阵 A 矩阵为 A = 抽样结果 T r μ = 服从二维正态分布 ( μ μ ),, σ σ σ T M 0 σ σ σ σ = σ σ σ σ η = μ + η σ = μ + [ + ] / σ η ( σ σ σ ) η y σ η y y 67