how to calculate similarity between any two points using custom distance fuction

I have a data set has m rows (data point) and 4 columns (properties) as below:

                i               j               k              p
1              30           50           30           1.2
2              30           60           30           1.3
…             ..             ..             ..             ..
m            80           40           40           0.2

i, j, k represent coordinates information, p represents other property say weight.

I want to define the distance or dissimilarity between any two data points as ( sqrt( (i1-i2)^2+(j1-j2)^2+(k1-k2)^2 ) / Dmax ) × (p1-p2 / Pmax), Where Dmax is the maximum distance between two points, and Pmax is the maximum weight gap between two points.

Can anyone give me a little tint for this one? If pdist is not working for this one, is there any other function that I can use? Or I have to write some code to calculate the dissimilarity every time, merge the points with smallest dissimilarity, update the dissimilarity matrix and original data matrix, merge, and do the circle. if this is the way, any efficient way to get it done?

Thanks a lot.

Derrick

[FirefoxMetzger]

What distance metric is Dmax based on?

Answer 1

[Star Strider]

I would leave out ‘Dmax’ and ‘Pmax’ because they’re constants and are not going to affect the relative distance results, and that is what you want. Those extra computations will simply impair the efficiency of your code.

I would create a distance function and use pdist:

M = [30           50           30           1.2
     30           60           30           1.3
     80           40           40           0.2];

dist_fcn = @(x1,x2) sqrt( (x1(1)-x2(1)).^2+(x1(2)-x2(2)).^2+(x1(3)-x2(3)).^2 ) .* abs(x1(4)-x2(4));

D = pdist(M, dist_fcn);
Dm = squareform(D)

Dm =

              0       1646.1       1646.1
         1646.1            0       60.249
         1646.1       60.249            0

See the documentation for pdist and squareform for details.

[Derrick Fu]

Thanks for your reply Star Strider. That dist_fcn makes sense to me, but if calculate distance between point 1 and 2: sqrt((30-30)^2+(50-60)^2+(30-30)^2).*|1.2-1.3| =0.316...<> 1646.1; do you have an idea about this or am I wrong somewhere? Thanks. Derrick

[Star Strider]

It took me a while to figure out what pdist wants its custom function to do. I did a hand calculation on the distance metrics, and this gives the correct answers:

M = [30           50           30           1.2
     30           60           30           1.3
     80           40           40           0.2];

dist_fcn = @(XI,XJ) sqrt(bsxfun(@minus,XI(:,1:3),XJ(:,1:3)).^2 * ones(3,1)) .* abs(XI(:,4)-XJ(:,4));
D = pdist(M, dist_fcn);

Ds = squareform(D)

Ds =

              0            1       51.962
              1            0       60.249
         51.962       60.249            0

By the way, I get:

sqrt((30-30)^2+(50-60)^2+(30-30)^2) * abs(1.2-1.3) = 1

since sqrt(10^2)*0.1 = 1.

So, at long last, it works. I learned much about how pdist works, so thank you for that opportunity.

[Derrick Fu]

Thanks a lot Star Strider. That is exactly what I am after. I compare the pdist and method mentioned by John BG using same custom distance calculation fuction, pdist is much quicker than the For...end way. Sorry for the late reply, but many thanks for the help.

[Derrick Fu]

Hi Star Strider, I have one more question if you could help. The "abs(XI(:,4)-XJ(:,4))" in the dist_fcn function you wrote, can we achieve something like this: IF(abs(XI(:,4)-XJ(:,4))=0, 0.001, abs(XI(:,4)-XJ(:,4))) . so the whole function is like dist_fcn = @(XI,XJ) sqrt(bsxfun(@minus,XI(:,1:3),XJ(:,1:3)).^2 * ones(3,1)) .* IF(abs(XI(:,4)-XJ(:,4))=0, 0.001, abs(XI(:,4)-XJ(:,4))); What I want to do is if abs(XI(:,4)-XJ(:,4))=0, it will be replaced by a very small value say 0.01, otherwise keep the result. Do you have an idea on this? any thoughts are welcome, I am quite new to matlab. Thanks in advance.

[Star Strider]

My pleasure.

The pdist and bsxfun functions are optimised, so they are both very fast.

I would use the max function:

dist_fcn = @(XI,XJ) sqrt(bsxfun(@minus,XI(:,1:3),XJ(:,1:3)).^2 * ones(3,1)) .* max(abs(XI(:,4)-XJ(:,4)), 0.01);

That runs (I tested it with the ‘M’ matrix here), and it produces the correct result. It takes the maximum value of ‘abs(XI-XJ)’ and 0.01, so it will do what you want.

M = [30           50           30           1.2
     30           60           30           1.3
     80           40           40           0.2
     20           40           10           0.2];

[Derrick Fu]

Oh my god, you are a genius! Thank you so much Star!!!

[Star Strider]

As always, my pleasure!

Answer 2

[John BG]

Hi Derrik

Your approach

( sqrt( (i1-i2)^2+(j1-j2)^2+(k1-k2)^2 ) / Dmax ) × (|p1-p2| / Pmax)

has a problem: while you input are, for each point, 4 variables

   [i j k p]

dividing by Dmax and Pmax is unnecessarily reducing resolution. So

1. arrange data in an m*4 array

ore1=[30           50           30           1.2;
        30           60           30           1.3;
         ..             ..             ..             ..
        80           40           40           0.2 ]

2. let's say you have m pairs to measure, arrange the pairs to measure in another array m*2

check1=[1 3;
               20 106;
               50 19;
                ..]
where each line of check(k,:) points at 2 different lines of ore1.

3. Since

 * the distance coordinates you are using do not seem to have negative values, 

 * and you are using the weight as a 4th dimension

try without measuring the euclidean distance, just with

     meas1=zeros(1,m)   % meas1 contains the measurements
     for k=1:1:m
             L1=orde1(check1(k,1),:);
             L2=orde1(check1(k,2),:);
             meas1(k)=sum(abs(L1-L2))
     end

If it's not the case that all coordinates are positive, perhaps you can shift the origin to a point far enough so that all coordinates are positive.

4. now in meas1 you have the pairs with higher similarity or closer each other in the ore with lower values, invert to have higher values proportional to higher similarity:

    meas12=max(meas1)-meas1;

cannot go much simpler than these lines.

Check the time with tic;toc; on a small sample to forecast how much time you may need to measure larger ores.

You can also try

(a*(i1-i2)^2+b*(j1-j2)^2+c*(k1-k2)^2+d*(p1-p2)^2)^.5

where a,b,c,d may be chosen as for instance

a=1;b=1;c=1;d=10

Try measuring the distance between 2 same points with different distance methods and you will realise that Dmax and Pmax are not needed.

You can still normalise the distance measurement, but after calculating distances.

If you find this answer of any help solving your question,

please click mark it as accepted answer on click on the thumbs-up vote link,

thanks in advance

John

[e-mail]

[Derrick Fu]

Thank you so much for the help really John.I appreciate that so much. Sorry I accept Star's answer, I think it would be helpful for others who have similar problems, please understand. Thanks a lot.

[John BG]

It's ok Derrik, you are the expert in data mining regards

John BG

<mailto:[e-mail] [e-mail]>

[Image Analyst]

Well, you can "Vote" for his answer (like I did) to throw a couple of reputation points his way to thank him for his efforts.

[Dina Abd El-twab]

If i have two circles , the red circle drawn using viscircles command , the blue circle drawn using (center(x,y) , raduis) of the circle . Then i want to calculate the dice similarity ( intersection ratio) between these two circles (red and blue) ,could you please help me ?

a=imread('US0018_0101.png');
figure, imshow(a)
%Guassian filtering 5x5
Iblur = imgaussfilt(a, 5);
 figure , imshow(Iblur)
 %morphological opening with disk of radius 9 
 sel=strel('disk',9);
d=imerode(Iblur,sel);
figure, imshow(d)
%thresholding with level of 0.15
level = 0.15;
b =im2bw(d,level);
figure ,imshow(b)
%hough transform to find circle(red)
Rmin = 25;
Rmax = 65;
[centersDark, radiiDark] = imfindcircles(b,[Rmin Rmax],'ObjectPolarity','dark');
viscircles(centersDark, radiiDark,'LineStyle','-');
hold on
% Draw blue circle with (x,y,r)
 h = circle(106,108,42)

@Image Analyst

Image Analyst

[Dina Abd El-twab]

[Dina Abd El-twab]

@Image Analyst

Image Analyst