title: Broadcast function in Numpy similar to matrix multiplication
body: <p>So I am building a KNN for some larger datasets, and I am required to run Leave-one-out cross validation in order to choose the correct K, so speed is important.  </p>
  
  <p>I am trying to do the distance calculation via broadcasting.  </p>
  
  <hr>
  
  <p>The situation is: <code>X</code> is my training matrix, a 2D matrix with samples on the rows.  <code>Q</code> is my query matrix, or test data, also with samples on the rows.  </p>
  
  <p>I need to run something analogous to a matrix multiplication, where I match each row of Q with each column of <code>X.T</code> (x transpose), and build a <code>sample x sample</code> matrix, where each entry [i,j] is the distance of query sample i from training sample j.  Then I would sort and choose the mode of the classes from the first k samples. </p>
  
  <hr>
  
  <p>Anyways, in matrix multiplication, numpy does exactly this...but instead of a distance calculation, it does a piecewise multiplication and sum (dot product).  If I could insert my distance function into that position instead, I think I would have a KNN distance calculation about as fast as a numpy matrix multiplication.  </p>
  
  <p>Is there a way to do this using broadcasting or some other numpy technique?</p>
  
  <p>Maybe even a way to parallelize it?</p>
  
  <hr>
  
  <p>Sample code: </p>
  
  <pre><code>import numpy as np
  
  x1 = np.asarray([1.0,10.0,100.0])
  x2 = np.asarray([40.0,60.0,80.0])
  x3 = np.asarray([20.,30.,40.])
  x = np.concatenate((x1.reshape(3,1),x2.reshape(3,1),x3.reshape(3,1)),axis=1)
  
  y1 = np.asarray([4.0,88.0,35.0])
  y2 = np.asarray([7.0,65.0,99.0])
  y3 = np.asarray([40.0,13.0,27.0])
  y = np.concatenate((y1.reshape(3,1),y2.reshape(3,1),y3.reshape(3,1)),axis=1)
  
  def euclidean_distance(p1,p2): 
      return np.sqrt(np.sum((p1-p2)**2.0))
  </code></pre>
  
  <p>so, I might write: </p>
  
  <pre><code>distances = np.zeros((y.shape[0],x.shape[0]))
  for i in range(y.shape[0]):
     for j in range(x.shape[0]):
         distances[i,j] = euclidean_distance(y[i,:],x[j,:])
  </code></pre>
  
  <p>And that is what I would then sort.  In the current for loop above, I would just choose my k nearest neighbors and find the class in that inner loop...but it would be much slower than calculating all the distances in a vectorized calculation.  </p>