Numpy/Scipy에서 평면으로부터 벡터의 직교 거리를 구하는 방법은 무엇입니까?

Question

나는 numpy 배열로 벡터 세트를 가지고있다. 두 벡터 v1과 v2에 의해 정의 된 평면으로부터 각각의 직교 거리를 구해야합니다. Gram-Schmidt 프로세스를 사용하여 단일 벡터에 대해이 작업을 쉽게 수행 할 수 있습니다. 각각의 벡터를 반복하거나 np.vectorize를 사용하지 않고도 많은 벡터에서 매우 빠르게 수행 할 수 있습니까?

감사합니다. 이 쉽게 할 수있는 세 가지 차원에서

:

Answer 1

제이미의 대답 @ 달성하기 위해 더 명시 적 방법, 당신은 프로젝션 연산자의 건설에 대한 명시 될 수있다 :

이

def normalized(v): 
    return v/np.sqrt(np.dot(v,v)) 
def ortho_proj_vec(v, uhat): 
    '''Returns the projection of the vector v on the subspace 
    orthogonal to uhat (which must be a unit vector) by subtracting off 
    the appropriate multiple of uhat. 
    i.e. dot(retval, uhat)==0 
    ''' 
    return v-np.dot(v,uhat)*uhat 

def ortho_proj_array(Varray, uhat): 
    ''' Compute the orhogonal projection for an entire array of vectors. 
    @arg Varray: is an array of vectors, each row is one vector 
      (i.e. Varray.shape[1]==len(uhat)). 
    @arg uhat: a unit vector 
    @retval : an array (same shape as Varray), where each vector 
       has had the component parallel to uhat removed. 
       postcondition: np.dot(retval[i,:], uhat) ==0.0 
       for all i. 
    ''' 
    return Varray-np.outer(np.dot(Varray, uhat), uhat) 




# We need to ensure that the vectors defining the subspace are unit norm 
v1hat=normalized(v1) 

# now to deal with v2, we need to project it into the subspace 
# orhogonal to v1, and normalize it 
v2hat=normalized(ortho_proj(v2, v1hat)) 
# TODO: if np.dot(normalized(v2), v1hat)) is close to 1.0, we probably 
# have an ill-conditioned system (the vectors are close to parallel) 



# Act on each of your data vectors with the projection matrix, 
# take the norm of the resulting vector. 
result=np.zeros(M.shape[0], dtype=M.dtype) 
for idx in xrange(M.shape[0]): 
    tmp=ortho_proj_vec(ortho_proj_vec(M[idx,:], v1hat), v2hat)    
    result[idx]=np.sqrt(np.dot(tmp,tmp)) 

# The preceeding loop could be avoided via 
#tmp=orhto_proj_array(ortho_proj_array(M, v1hat), v2hat) 
#result=np.sum(tmp**2, axis=-1) 
# but this results in many copies of matrices that are the same 
# size as M, so, initially, I prefer the loop approach just on 
# a memory usage basis. 

이 정말 그램 - 슈미트 직교 절차의 단지 일반화이다. 참고이 과정의 끝에서 우리가 np.dot(v2hat,v2hat.T)==1, np.dot(v1hat, v2hat)==0

Answer 2

당신은 단위 정상 평면을 구성 할 필요가

nhat=np.cross(v1, v2) 
nhat=nhat/np.sqrt(np.dot(nhat,nhat)) 

다음 벡터의 각각이 점; 이는 I 가정이되도록 result[idx] 비행기에서 idx'th 벡터의 거리 인 Nx3 매트릭스 M

result=np.zeros(M.shape[0], dtype=M.dtype) 
for idx in xrange(M.shape[0]): 
    result[idx]=np.abs(np.dot(nhat, M[idx,:])) 

이다.

Answer 3

EDIT 내가 작성한 원본 코드가 올바르게 작동하지 않아 제거했습니다.

def distance(v1, v2, u) : 
    u = np.array(u, ndmin=2) 
    v = np.vstack((v1, v2)) 
    vv = np.dot(v, v.T) # shape (2, 2) 
    uv = np.dot(u, v.T) # shape (n ,2) 
    ab = np.dot(np.linalg.inv(vv), uv.T) # shape(2, n) 
    w = u - np.dot(ab.T, v) 
    return np.sqrt(np.sum(w**2, axis=1)) # shape (n,) 

가 작동하는지 확인하려면 다음을하지만 같은 생각 다음, 당신이 그것을 생각하는 시간을 보낼 경우, 크레이머의 규칙에 대한 필요가 없습니다, 다음과 같이 코드가 상당히 간소화 될 수 있으며, 아래에 설명 제대로, 나는 distance_3d로 함수에 데이브의 코드를 포장하고 다음을 시도했다 :

>>> d, n = 3, 1000 
>>> v1, v2, u = np.random.rand(d), np.random.rand(d), np.random.rand(n, d) 
>>> np.testing.assert_almost_equal(distance_3d(v1, v2, u), distance(v1, v2, u)) 

그러나 물론 지금은 어떤 d 작동 :

>>> d, n = 1000, 3 
>>> v1, v2, u = np.random.rand(d), np.random.rand(d), np.random.rand(n, d) 
>>> distance(v1, v2, u) 
array([ 10.57891286, 10.89765779, 10.75935644]) 

w 따라서 np.dot(w, v1) = np.dot(w, v2) = 0면에 수직이고, 반면 94,922,398,581,088,243,210 당신은 당신의 벡터를 분해해야

는 v 비행기에, 그래서 v = a * v1 + b * v2로 분해 될 수있다, 두 벡터, u = v + w의 합에, u를 호출 할 수 있습니다 . 당신이 u = a * v1 + b * v2 + w를 작성하고 v1 및 v2이 표현의 내적을 취하면

, 당신은 두 개의 미지수 두 개의 방정식을 얻을 : 그것은 단지 2 × 2 시스템이기 때문에

np.dot(u, v1) = a * np.dot(v1, v1) + b * np.dot(v2, v1) 
np.dot(u, v2) = a * np.dot(v1, v2) + b * np.dot(v2, v2) 

을, 우리는 Cramer's rule를 사용하여 해결할 수 있습니다 로 : 여기에서

uv1 = np.dot(u, v1) 
uv2 = np.dot(u, v2) 
v11 = np.dot(v1, v2) 
v22 = np.dot(v2, v2) 
v12 = np.dot(v1, v2) 
v21 = np.dot(v2, v1) 
det = v11 * v22 - v21 * v12 
a = (uv1 * v22 - v21 * uv2)/det 
b = (v11 * uv2 - uv1 * v12)/det 

, 당신은 얻을 수 있습니다 :

w = u - v = u - a * v1 - b * v2 

이고 평면까지의 거리는 모듈 w입니다.