我使用python与肥胖型(约2000×2000)的矩阵,其中每个工作I,J在矩阵点代表单个像素.
矩阵本身是稀疏的(即它们的大部分将具有零值),但是当它们被更新时,它们倾向于增量操作,对于矩形"块"中的大量相邻像素,而不是这里的随机像素或那里(我目前没有利用的财产......).对矩阵算法有点新意,但我已经研究了许多可能的解决方案,包括scipy稀疏矩阵的各种风格.到目前为止,协调(COO)矩阵似乎是最有希望的.因此,例如,我想增加一个块形状,我必须做的事情是:
>>> from scipy import sparse
>>> from numpy import array
>>> I = array([0,0,0,0])
>>> J = array([0,1,2,3])
>>> V = array([1,1,1,1])
>>> incr_matrix = sparse.coo_matrix((V,(I,J)),shape=(100,100))
>>> main_matrix += incr_matrix #where main_matrix was previously defined
在未来,我希望在任何情况下都有更丰富的像素值表示(代表RGB等的元组),numpy数组不支持开箱即用(或者我可能需要使用它).
最终我会有一些这些矩阵,我需要做简单的算术,我需要代码尽可能高效 - 并且可以分发,所以我需要能够坚持和交换这些物体在小型表示中没有实质性的惩罚.我想知道这是否是正确的方法,或者我是否应该使用dicts等自己的结构？
我正在尝试使用scipy的稀疏包来实现以下等式:
W = x[:,1] * y[:,1].T + x[:,2] * y[:,2].T + ...
其中x和y是nxm csc_matrix.基本上我试图将x的每个col乘以y的每个col并将得到的nxn矩阵求和.然后我想制作所有非零元素1.
这是我目前的实施:
c = sparse.csc_matrix((n, n))
for i in xrange(0,m):
tmp = bam.id2sym_thal[:,i] * bam.id2sym_cort[:,i].T
minimum(tmp.data,ones_like(tmp.data),tmp.data)
maximum(tmp.data,ones_like(tmp.data),tmp.data)
c = c + tmp
此实现具有以下问题:
内存使用似乎爆炸了.根据我的理解,内存应该只会随着c变稀疏而增加,但我看到循环开始吞噬> 20GB的内存,其中= 10,000,m = 100,000(x和y的每一行只有大约60个非零元素).
我正在使用一个效率不高的python循环.
我的问题:有更好的方法吗？控制内存使用量是我首先关注的问题,但要让它更快更好!
谢谢!
假设我有一个非常大的稀疏数据矩阵,但我只想看一下它的样本,使它更稀疏.假设我还有一个三元数据框,包括数据的行/列/值列(从csv文件导入).我知道我可以使用库(Matrix)的sparseMatrix()函数来创建稀疏矩阵
sparseMatrix(i=df$row,j=df$column,x=df$value)
但是,由于我的值,我最终得到了一个稀疏矩阵,其中包含数百万行数万列(其中大多数都是空的,因为我的子集排除了大多数行和列).所有这些零行和列最终都会扭曲我的一些函数(例如,采用聚类 - 当一个原点甚至不是有效点时,我最终会得到一个包含原点的聚类).我想执行相同的操作,但使用i和j作为rownames和colnames.我已经尝试创建一个密集的矢量,采样到最大尺寸并使用添加值
denseMatrix <- matrix(0,nrows,ncols,dimnames=c(df$row,df$column))
denseMatrix[as.character(df$row),as.character(df$column)]=df$value
(实际上我一直把它设置为等于1,因为我对这种情况下的值不感兴趣)但是我发现它填充了整个矩阵,因为它需要所有行和列的交叉而不仅仅是row1*col1,row2*col2 ...有没有人知道如何完成我想要做的事情？或者我可以填写一个稀疏矩阵并简单地以某种方式丢弃所有零行和列以将其自身压缩为更密集的形式(但我想保留一些参考回原始行和列号)我感谢任何建议!
这是一个例子:
> rows<-c(3,1,3,5)
> cols<-c(2,4,6,6)
> mtx<-sparseMatrix(i=rows,j=cols,x=1)
> mtx
5 x 6 sparse Matrix of class "dgCMatrix"
[1,] . . . 1 . .
[2,] . . . . . .
[3,] . 1 . . . 1
[4,] . . . . . .
[5,] . . . . . 1
我想摆脱colums 1,3和5以及第2行和第4行.这是一个非常简单的例子,但想象一下,如果没有行号1,3和5,它们分别是1000,3000和5000.然后它们之间会有更多的空行.这是当我使用具有命名行/列的密集矩阵时发生的情况
> dmtx<-matrix(0,3,3,dimnames=list(c(1,3,5),c(2,4,6)))
> dmtx
2 4 6
1 0 0 0
3 0 0 0 …我试图在Mathematica中解决一个相当大的线性规划问题,但由于某种原因,瓶颈是设置线性约束的数组.
我初始化矩阵的代码如下所示:
AbsoluteTiming[S = SparseArray[{{i_, 1} -> iaa[[i]],
{i_, j_} /; Divisible[a[[j]], aa[[i]]] -> 1.}, {2*n, n}]]
这里n是4455,iaa是实数列表,a,aa是整数列表.我得到的这条线的输出是
{2652.014773,SparseArray[<111742>,{8910,4455}]}
换句话说,初始化此矩阵需要45分钟,即使它只有111,742个非零项.为了比较,实际上解决线性程序只需要17秒.是什么赋予了？
编辑:此外,任何人都可以解释为什么这会占用如此多的内存,因为它运行？因为在用户时间,这个计算花费的时间不到十分钟......大部分计算时间都花费在内存中进行分页.
Mathematica出于某种原因将这个矩阵存储为非稀疏矩阵,而它正在构建它吗？因为那真的很蠢.
好吧,我不认为,我可以用文字解释这个问题所以,这里是ipython会话的片段,我导入scipy,以构建一个稀疏矩阵.
In [1]: import scipy as sp
In [2]: a = sp.sparse.lil_matrix((5,5))
---------------------------------------------------------------------------
AttributeError Traceback (most recent call last)
/home/liveuser/<ipython-input-2-b5a55fc2d0ac> in <module>()
----> 1 a = sp.sparse.lil_matrix((5,5))
AttributeError: 'module' object has no attribute 'sparse'
In [3]: import scipy.sparse as spar
In [4]: ax = spar.lil_matrix((5,5))
In [5]: a = sp.sparse.lil_matrix((5,5)) # you are kidding me?
In [6]: a
Out[6]:
<5x5 sparse matrix of type '<type 'numpy.float64'>'
with 0 stored elements in LInked List format>
In [7]: ax
Out[7]:
<5x5 sparse …我有一个稀疏的数组,似乎太大了,无法在内存中有效地处理handel(2000x2500000,float).我可以把它形成一个稀疏的lil_array(scipy),但如果我尝试输出一个列或行压缩的稀疏数组(A.tocsc(),A.tocsr())我的机器内存不足(并且之间也存在严重的不匹配)文本文件4.4G和pickeled lil数组12G中的数据 - 拥有一个更接近原始数据大小的磁盘格式会很不错.
我将来可能会处理更大的阵列.
问题:处理大型磁盘阵列的最佳方法是什么,以便我可以透明地使用常规numpy函数.例如,行和列的总和,矢量积,最大值,最小值,切片等？
是什么方式去pytables？是否有一个好的(快速)sql-numpy中间件层？内置于numpy的磁盘阵列上的秘密？
在过去使用(稍微小一点)的数组时,我总是将长期计算结果缓存到磁盘上.当阵列最终<4G左右但不再成立时,这种方法有效.
我有这个问题,这需要解决X在AX=B.A是15000 x 15000的顺序,是稀疏和对称的.B是15000 X 7500并且不稀疏.解决X的最快方法是什么？
我可以想到两种方式.
X = A\B使用for循环,
invA = A\speye(size(A))
for i = 1:size(B,2)
X(:,i) = invA*B(:,i);
end
有没有比上面两个更好的方法？如果没有,我提到的两者之间哪一个最好？
matlab linear-algebra sparse-matrix matrix-inverse factorization
我正在处理一个非常大的稀疏矩阵乘法(matmul)问题.举个例子,让我们说:
A是二进制(75 x 200,000)矩阵.它很稀疏,所以我使用csc进行存储.我需要做以下matmul操作:
B = A.transpose()*A.
输出将是一个大小为200Kx200K的稀疏对称矩阵.
不幸的是,B将是方式在我的笔记本电脑到大的RAM来存储(或"核心").另一方面,我很幸运,因为B有一些属性可以解决这个问题.
由于B将沿对角线和稀疏对称,我可以使用三角矩阵(上/下)来存储matmul操作的结果,稀疏矩阵存储格式可以进一步减小尺寸.
我的问题是......可以提前告诉numpy或scipy,输出存储要求会是什么样子,以便我可以使用numpy选择存储解决方案并避免"矩阵太大"运行时错误几分钟(小时)的计算？
换句话说,通过使用近似计数算法分析两个输入矩阵的内容,可以近似地对矩阵的存储要求进行近似.
如果没有,我正在寻找一个强力解决方案.来自以下Web链接的涉及map/reduce,核心外存储或matmul细分解决方案(strassens算法)的内容:
一对Map/Reduce问题细分解决方案
一种核外(PyTables)存储解决方案
一个matmul细分解决方案:
提前感谢任何建议,评论或指导!
python numpy combinatorics sparse-matrix matrix-multiplication
我有一个大型矩阵存储为scipy.sparse.csc_matrix,并希望从大矩阵中的每一列中减去一个列向量.当你正在进行规范化/标准化这样的事情时,这是一个非常常见的任务,但我似乎找不到有效地做到这一点的正确方法.
这是一个示例:
# mat is a 3x3 matrix
mat = scipy.sparse.csc_matrix([[1, 2, 3],
[2, 3, 4],
[3, 4, 5]])
#vec is a 3x1 matrix (or a column vector)
vec = scipy.sparse.csc_matrix([1,2,3]).T
"""
I want to subtract `vec` from each of the columns in `mat` yielding...
[[0, 1, 2],
[0, 1, 2],
[0, 1, 2]]
"""
实现我想要的一种方法是对vec自己进行3次自身,产生一个3x3矩阵,其中每列都是vec,然后从中减去mat.但同样,我正在寻找一种有效地做到这一点的方法,并且hstacked矩阵需要很长时间才能创建.我确信有一些神奇的方法可以用切片和广播来做到这一点,但它让我望而却步.
谢谢!
编辑:删除了"就地"约束,因为稀疏结构将在就地分配方案中不断变化.
我一直在使用R中的库'doParallel'来提高一组函数的速度.但是,我遇到了一个我无法解决的错误.我相信以下代码隔离了问题的精髓:
library(Matrix)
library(doParallel)
test_mat = Matrix(c(0,1,2,NA,0,0,2,NA,1,NA,1,2,2,NA,0,1,0,2,2,2,0,0,NA,NA,1,2,1,1,2,1,rep(NA,5)), ncol=7, byrow=TRUE, sparse=TRUE)
par_func <- function(mat, ncores)
{
cl <- makePSOCKcluster(ncores)
clusterSetRNGStream(cl)
registerDoParallel(cl, cores = ncores)
df = data.frame(1:7, NA)
temp_vec = foreach(i=iter(df, by='row'), .combine=rbind) %dopar%
{
i[,2] <- sum(mat[,i[,1]] == 1, na.rm = TRUE) + 1
}
stopCluster(cl)
return(temp_vec)
}
par_func(mat=test_mat, ncores=5)
这会产生以下错误消息:
Error in { : task 1 failed - "object of type 'S4' is not subsettable"
如果'mat'是'matrix'类而不是'dgCMatrix',则此函数有效,因此问题似乎是由于稀疏矩阵的子集化.我有什么选择可以解决这个问题吗？矩阵"mat"可以非常大并且可以包含许多零,因此我想继续使用稀疏矩阵.