作者alen332l (alen3321)
看板Python
标题Re: [问题] 关於反矩阵inv计算爆值是正常的吗?
时间Thu Jun 23 11:02:36 2016
※ 引述《alen332l (alen3321)》之铭言：
: ※ 引述《ray0215 (❺❺❻❻)》之铭言：
: : 小弟我最近需要计算相关系数(correlation coefficient matri)矩阵的反矩阵
: : 例如
: : Martrix array([[ 1. , 0.92531576, 0.81998247, 0.92586894],
: : [ 0.92531576, 1. , 0.70974912, 0.72358554],
: : [ 0.81998247, 0.70974912, 1. , 0.72669914],
: : [ 0.92586894, 0.72358554, 0.72669914, 1. ]])
: : 这样用numpy.linalg.inv算出来的值
: : 会计算出
: : array([[ 6.51500296e+15, -3.18743561e+15, -7.89299760e+14,
: : -3.15207310e+15],
: : [ -3.18743561e+15, 1.55943840e+15, 3.86161323e+14,
: : 1.54213745e+15],
: : [ -7.89299760e+14, 3.86161323e+14, 9.56245323e+13,
: : 3.81877116e+14],
: : [ -3.15207310e+15, 1.54213745e+15, 3.81877116e+14,
: : 1.52502845e+15]])
: : 每个都十五次方
: : 我想应该不是正确解答吧@@
: : 请问要怎麽解决
: : 先谢谢版上的大大了
: 你推文说要计算R^(-1)-I 的值
: R = correlation matrix
: I = idendity matrix
: 处理矩阵时
: 有一个原则
: 能不要反矩阵 就不要算 因为算反矩阵很昂贵 又不稳定
: 你原本题目应该是：
: R ( X + I ) = I ,其中X为你要求的矩阵.
: 那麽，这样处理比较稳定:
: R X = (I - R)
: └已知┘ └已知┘
: X 切成column表达:
: X = [x1 │x2 │x3...] ,其中x1为向量,x2,x3...亦同
: (I-R)定义成 B 好了
: 也切成column:
: (I-R) =def= B = [b1│b2│b3...]
: 原问题变成
: R ‧ x1 = b1
: (nxn) (nx1) (nx1) →维度验算
: R ‧ x2 = b2
: R ‧ x3 = b3
: 每个问题用 BLAS 算完後
: 再把 x1,x2,x3...拼回来
: 印象中是这样,几年前处理的印象...
: p.s.
: 1.回文如有错请高手指正 已经有一段时间没处理matrix了 谢谢
: 2.我没实际算...如果算出来是合理值 那就是你原本用inv时会有误差造成
: 行列式小时，算inv的误差会累积
: 通常基础线性代数不会讲
: 数学系的纯数课程也不会讲
我把我刚才的回文写成code
如下: 但答案还是有10的7次方就是了
----------------------------------------------------------------
import numpy as np
R =np.array([[ 1. , 0.92531576, 0.81998247, 0.92586894],
[ 0.92531576, 1. , 0.70974912, 0.72358554],
[ 0.81998247, 0.70974912, 1. , 0.72669914],
[ 0.92586894, 0.72358554, 0.72669914, 1. ]])
I = np.array([[1,0,0,0],[0,1,0,0],[0,0,1,0],[0,0,0,1]])
B = I - R
print("R",R)
print("B",B)
X = np.linalg.solve(R, B)
print("X",X)
-----------------------------------------------
答案就是 X
--
※ 发信站: 批踢踢实业坊(ptt.cc), 来自: 36.228.190.212
※ 文章网址: https://webptt.com/cn.aspx?n=bbs/Python/M.1466650959.A.B2B.html
※ 编辑: alen332l (36.228.190.212), 06/23/2016 18:47:27