I'm having the same problem

Python 3.6.6 |Anaconda, Inc.| (default, Jun 28 2018, 11:27:44) [MSC v.1900 64 bit (AMD64)] on win32

Same NumPy, SciPy, sklearn, Windows 10 build 17134.

Thanks for the report!

See #9354 There is a know issue with the precision of pairwise_distances when using 32 bit floats that we really want to fix. A solution was proposed in #11271 but it's still under discussion. If any of you have thoughts on the linked issue and Pull request you are very welcome to comment there.

I think it's due to the above linked issue, but I also can't reproduce the results of your code snippet with a comparable configuration on Linux running it multiple times with 0.19.1 and master. If you could define your random datasets with a RNG, (e.g. np.random.RandomState(42).rand(9, 10000).astype('float32')) that might help.

I modified the snippet to use a seed.

Note that when using float64, the error is -7.815970093361102e-14, so this is not exactly the same as #9354. The problem is more strange than #9354 because it means the computation with the same function doesn't give twice the same result.

The context of this code is that I compute the distance matrix of a high-dimensional dataset, so I proceed by chunks and fill small squares of the matrix. I noticed that the result depends on the chunk size, which is not normal.

@jnothman I'm not complaining about a numerical precision issue. I'm complaining about the fact that the result depends on the other elements of the matrix.

I see. Sorry, haven't had my coffee. Interesting. Deserves investigation.

I can reproduce. It seems to come from NumPy:

import numpy as np
np.random.seed(42)

a = np.array(np.random.rand(9, 10000), dtype=np.float32)
m1 = np.dot(a, a.T)
m2 = np.dot(a[:4], a[:4].T)

print(m2[0, 1] - m1[0, 1])
# 0.0012207

Linux-3.16.0-6-amd64-x86_64-with-debian-8.11
Python 3.6.2 |Continuum Analytics, Inc.| (default, Jul 20 2017, 13:51:32) 
[GCC 4.4.7 20120313 (Red Hat 4.4.7-1)]
NumPy 1.13.1
SciPy 0.19.1
Scikit-Learn 0.20.dev0

Note that this can probably be considered as a normal numerical error (see e.g. numpy/numpy#11419 (comment)).
Indeed, np.finfo(np.float32).eps * 10000 = 0.0011920 is of the same order as the error 0.0012207.

Thanks for the code snipped @TomDLT ! Reported the issue upstream in numpy/numpy#11655 with additional details.

It might be linked to MKL, might be also the expected ~7 significant digits precision of 32 bit floats. In any case, this is a fairly low level issue and numpy would be a better place to discuss it.

Should we close this in favor of the numpy issue?

I can reproduce the snippet of @TomDLT, but it prints -0.00048828125.

What is strange is that some machines don't have this problem.

I agree to close this here!

Thanks for the detailed report in any case @louisabraham !