siavashk / pycpd

Pure Numpy Implementation of the Coherent Point Drift Algorithm
MIT License
513 stars 115 forks source link

Fixed typo in em_registration _init_ and error in gaussian_kernel() #23

Closed gattia closed 5 years ago

gattia commented 5 years ago

Assuming there wasnt an error in the original CPD paper (I looked at two versions of it) it seems there is a type in the guassian_kernel(), the beta coefficient should be squared.

And then in the deformable_registration() class __init__() function it seems that there might have been a typo where beta is set to 2 (default) based on if alpha is None. It seems this should have been based on if beta is None. I might be missing something here though because after noticing this, it also seemed to me like setting self.alpha and self.beta in this way was unnecessary as they are already set as defaults in the deformable_registration() function. So, if they arent prescribed by the user calling the function they should already be set to the default of 2. So having them at all (lines 19 and 20 in deformable_registration.py) might be redundant.

gattia commented 5 years ago

I've mixed something up in here... sorry. I was only attempting to merge the first two commits:

[ea65000](/siavashk/pycpd/pull/23/commits/ea65000b49eb4b27a2af3662cc08dc2f65d8d3e2) [40ba9ff](/siavashk/pycpd/pull/23/commits/40ba9fff14f933fa8caec37a33db5e5de104c809)

The others somehow got in the mix. I'll try to clean up.

Edit. If some of the other changes (mostly simplifying some of the functions) are of interest, then they were something I was going to propose once I'd played around with them a bit more on my end.

gattia commented 5 years ago

Ignore my prior comment, I sorted it out and got rid of the extra commits.

If there is interest in some other syntactical changes (that I mentioned in my prior comment), I'd be happy to do another pull request.

For example,

def initialize_sigma2(X,Y):
    (N, D) = X.shape
    (M, _) = Y.shape
    XX = np.reshape(X, (1,N,D))
    YY = np.reshape(Y, (M,1,D))
    XX = np.tile(XX, (M,1,1))
    YY = np.tile(YY, (1,N,1))
    diff = XX-YY
    err = np.multiply(diff, diff)
    return np.sum(err) / (D*M*N)

can be simplified to:

def initialize_sigma2(X,Y):
    (N, D) = X.shape
    (M, _) = Y.shape
    diff = X[None,:,:] - Y[:,None,:]
    err = diff**2
    return np.sum(err) / (D*M*N)

This same use of None indexing also works to remove the for loop that is currently being implemented in expectations(). The current loop could be replaced with:

P = np.sum((self.X[None,:,:] - self.TY[:,None,:])**2, axis=2)

siavashk commented 5 years ago

Thank you for your contributing to this repository. I just got back from vacation last night. I will go through your commits in the next couple of days.

siavashk commented 5 years ago

If there is interest in some other syntactical changes (that I mentioned in my prior comment), I'd be happy to do another pull request.

Please do.