1

I am reading the same binary file into Python and Matlab and putting it into a matrix. When I take the norm of this matrix, I get different results.

I am using the respective smatload function to load the binary file.

Python:

def smatload(filename):
    #print 'opening: ', filename
    f = open(filename, 'rb')
    m = np.fromfile(f,'q',1)
    n = np.fromfile(f,'q',1)
    nnz = np.fromfile(f,'q',1)
    print 'reading %d x %d with %d non-zeros' % (m,n,nnz)
    S = np.fromfile(f,'d',3*nnz)
    f.close()
    S = S.reshape((nnz,3))
    rows = S[:,0].astype(int) - 1
    cols = S[:,1].astype(int) - 1
    vals = S[:,2]
    return csr_matrix((vals,(rows,cols)),shape=(m,n))

Matlab:

function [A] = smatload(filename)

fid = fopen(filename,'r');
if( fid == -1 )
    disp(sprintf('Error: Unable to open file [%s], fid=%d\n',filename,fid));
    A = [-1];
    fclose(fid);
    return;
end

m   = fread(fid,[1 1],'uint64');
n   = fread(fid,[1 1],'uint64');
nnz = fread(fid,[1 1],'uint64');

fprintf('Reading %d x %d with %d non-zeros\n',m,n,nnz);

S = fread(fid,[3 nnz],'double');
fclose(fid);
A = sparse(S(1,:)',S(2,:)',S(3,:)',m,n);

The results that I get for the norms of the returned matrices are

Matlab: norm(A.'fro') = 0.018317077159881

Python: np.linalg.norm(A) = 0.018317077159760

I have confirmed that they are both reading the correct number of values (6590x7126 matrix, 122526 non-zeros), and that I am using the same norm for both (frobenius).

Any ideas as to what might cause this?

Improve this question
12
  • This could just be a case of MATLAB doing various things to increase accuracy of floating point mathematics. I know that in MATLAB, if you add 0.01 1,000,000 times, you'll actually get 10,000, rather than the slightly off version you get in most other languages. Commented Feb 24, 2015 at 23:25
  • I have used this function other times for projects that I am doing and I haven't seen a difference. I'm curious as to why this is popping up now. It is causing an issue that my Matlab solution did not run into. Commented Feb 24, 2015 at 23:27
  • So you are concerned about differences off in the 12th decimal place when dealing with 100,000 floats? Not that it's relevant, but I'd be curious if the MATLAB array can be saved to a .mat file, and read with scipy.io.loadmat, and whether the results are the same. Commented Feb 24, 2015 at 23:35
  • Given the experiments that I am doing, it is causing some issues with convergence criterion. Let me check what that gives, and I'll update my post. Commented Feb 24, 2015 at 23:38
  • @David did your previous uses of the function involve such large matrices? I could easily see a 6590x7126 matrix getting into numbers too large for typical floating point representations to have high accuracy. Commented Feb 24, 2015 at 23:40

3 Answers 3

Reset to default
4

A quick look at the Frobenius Norm shows it requires squaring all of the values and adding them together.

Since you have uint64 in the read command, it looks like you could be filling up the floating point storage. When you multiply two binary numbers, it takes twice as many bits to store the answer. This means you would need 128 bits to store all the decimal values. If Python and MATLAB are doing this differently, that could explain why your decimal values are different.

See these two links for information about how MATLAB and Python handle floating point accuracy:

Python: https://docs.python.org/2/tutorial/floatingpoint.html

MATLAB: http://blogs.mathworks.com/cleve/2014/07/07/floating-point-numbers/

Improve this answer
Sign up to request clarification or add additional context in comments.

8 Comments

Not quite sure how to upload a file here, but posting the numbers isn't quite possible.
Also, have you confirmed that you are getting the same values going into the separate calculations?
Yes, I have. I use the norm of this matrix for some of the calculations, so this difference is propagating throughout the program.
Could you just add the numbers that produced the value '0.018317077159881' and '0.018317077159760'?
I'll try a sum, one second.
|
2

Well Matlab definitely seems to have different implementation for sparse and for dense arrays. Using the 4425x7126 sparse matrix A with 54882 non zero entries that you linked to and the following commands : 

FA=full(A);
av=A(:);
fav=FA(:);

I would expect the following commands to all yield the same value, because they are all computing the square root of the sum of the squares of the (non zero) elements of A : 

norm(A,'fro')
norm(av,2)
norm(FA,'fro')
norm(fav,2)

sqrt( sum(av .* av) )
sqrt( sum(av .^ 2) )

sqrt( sum(fav .* fav) )
sqrt( sum(fav .^ 2) )

In fact we see three slightly different answers : 

 norm(A,'fro')             0.0223294051001499
 norm(av,2)                0.0223294051001499
 norm(FA,'fro')            0.0223294051001499
 norm(fav,2)               0.0223294051001499

 sqrt( sum(av .* av) )     0.0223294051001521
 sqrt( sum(av .^ 2) )      0.0223294051001521

 sqrt( sum(fav .* fav) )   0.0223294051001506
 sqrt( sum(fav .^ 2) )     0.0223294051001506

In fact, even the reported sum of the elements of the sparse and dense representations of A are (a little bit) different : 

sum(A(:))                 1.00000000000068
sum(FA(:))                1.00000000000035

These differences seem to be of the same order of magnitude that you were seeing between the Python and Matlab norms.

Improve this answer

Comments

1

Not an answer but I don't have enough rep to comment. Would it be worthwhile to try to narrow the problem down a bit? If you partitioned your original matrix into 4 sub-matrices (Top left, top right, bottom left, bottom right) and compared the Frobenius norm reported in Matlab and in Python for each sub-matrix do you still see a discrepency between any values? If yes, then rinse and repeat on that sub-matix. If no, then don't even waste your time reading this comment. :-)

Improve this answer

9 Comments

I can try that. I'm worried that since the norm is off, but the original data matches, that it may be a difference in how matlab/numpy implement sqrt and ^2. If so then I won't be able to fix the issue.
I guess your original matrix must be about 375 MB on disk. If you could send it to me somehow (dropbox, etc) I could run your sample Matlab and Python code on my machine to calculate the norms. That might tell you if there is something version/environment specific causing the difference. (Of course it's more likely that we'd just have four different values for the norm - two different ones from your machine and two more different ones from mine. :-)
Are your entries real or complex? How about if you just compute the norm yourself? Or the squared norm? Do you calculations agree with the corresponding Matlab and Python versions? Does is make a difference if you use exponentiation vs plain old multiplication in calculating the (squared) norms? You mentioned that one of the implementations (I forget which) gave the "correct' answer - how do you know? (It's not that I am doubting you, I just don't know how you know)
I'm assuming that matlab is correct as it runs to completion and runs correctly. Python errors out halfway during an svd computation as i am using the norm to calculate how many values I want. I'll upload the files and send you a link in a minute.
The file is here The size is 1.3mb as everything is encoded in binary. The format is unsigned 64 bit int for the first dimension, unsigned 64 bit int for second dimension, unsigned 64 bit int for number of non-zeros, then number of non-zero entries consisting of three 64 bit doubles of row, col, value. I will try your suggestions but I did do the norm manually in python and got the same results. I even computed the square root manually to below machine precision.
|

Start asking to get answers

Find the answer to your question by asking.

Ask question

Explore related questions

See similar questions with these tags.