BookmarkSubscribeRSS Feed
🔒 This topic is solved and locked. Need further help from the community? Please sign in and ask a new question.
SWEETSAS
Obsidian | Level 7

 

Here are two data sets: 

Data have1;

Input obs reps id $ x  y;

Datalines ;

 

 

Obs

REPS

ID

X

Y

1

1

00-01

5

3

2

1

00-02

4

6

3

1

00-03

6

4

4

1

00-06

4

6

5

2

00-02

4

6

6

2

00-03

6

4

7

2

00-04

7

5

8

2

00-05

9

2

9

3

00-01

5

3

10

3

00-02

4

6

11

3

00-05

9

2

12

3

00-10

7

4

13

4

00-03

6

4

14

4

00-05

9

2

15

4

00-09

7

3

16

4

00-10

7

4

17

5

00-03

6

4

18

5

00-04

7

5

19

5

00-07

5

8

20

5

00-09

7

3

21

6

00-02

4

6

22

6

00-03

6

4

23

6

00-05

9

2

24

6

00-10

7

4

 

 

 

Data have2;

Input obs reps id $ x y;

Datalines ;

25

7

00-01

5

3

26

7

00-08

1

8

27

7

00-09

7

3

28

7

00-10

7

4

 

;

Run;

 

Dataset HAVE1 has 6 replicates indicated by variable REPS. I want to compute the Mahalanobis distance between the center of HAVE2 and each of the REPS in HAVE1. The vector of means of x and y in each of the datasets can serve as the center of each of the datasets. The distance will then added to columns in HAVE1. Does any know how to proceed?

 

Thanks in advance!

 

Jack

 

1 ACCEPTED SOLUTION

Accepted Solutions
Ksharp
Super User
OK. This code you don't need to consider about the number of variables.
Just define the calculated variables(x y ........ ) at :

read all var {x y .........} into data[c=vnames];

and Here :
read all var {x y ........} into m;





Data have1;
Input obs reps id $ x  y;
Datalines ;
1
1
00-01
5
3
2
1
00-02
4
6
3
1
00-03
6
4
4
1
00-06
4
6
5
2
00-02
4
6
6
2
00-03
6
4
7
2
00-04
7
5
8
2
00-05
9
2
9
3
00-01
5
3
10
3
00-02
4
6
11
3
00-05
9
2
12
3
00-10
7
4
13
4
00-03
6
4
14
4
00-05
9
2
15
4
00-09
7
3
16
4
00-10
7
4
17
5
00-03
6
4
18
5
00-04
7
5
19
5
00-07
5
8
20
5
00-09
7
3
21
6
00-02
4
6
22
6
00-03
6
4
23
6
00-05
9
2
24
6
00-10
7
4
;
run;
Data have2;
Input obs reps id $ x y;
Datalines ;
25
7
00-01
5
3
26
7
00-08
1
8
27
7
00-09
7
3
28
7
00-10
7
4
;
Run;

proc iml;
use have1 nobs nobs;
read all var {reps};
read all var {x y} into data[c=vnames];
close;
use have2;
read all var {x y} into m;
close;

center=m[:,];
level=unique(reps);
want=j(ncol(level),ncol(data)+2,.); 

do i=1 to ncol(level);
 want[i,1]=level[i];
 idx=loc(reps=level[i]);
 want[i,2:ncol(data)+1]=data[idx,][:,];
end;
 
 xx=want[,2:ncol(data)+1];
 cov=cov(xx);
 want[,ncol(data)+2]= mahalanobis(xx,center,cov);

names={reps}||vnames||{distance};
create want from want[c=names];
append from want;
close;

quit;

View solution in original post

23 REPLIES 23
Reeza
Super User

Break it into parts.

 

1. Calculate centroids - proc means

2. Merge datasets 

3. Calculate distance

Ksharp
Super User
So means of x and y in HAVE2 is the center of Mahalanobis distance ?
Data have1;
Input obs reps id $ x  y;
Datalines ;
1
1
00-01
5
3
2
1
00-02
4
6
3
1
00-03
6
4
4
1
00-06
4
6
5
2
00-02
4
6
6
2
00-03
6
4
7
2
00-04
7
5
8
2
00-05
9
2
9
3
00-01
5
3
10
3
00-02
4
6
11
3
00-05
9
2
12
3
00-10
7
4
13
4
00-03
6
4
14
4
00-05
9
2
15
4
00-09
7
3
16
4
00-10
7
4
17
5
00-03
6
4
18
5
00-04
7
5
19
5
00-07
5
8
20
5
00-09
7
3
21
6
00-02
4
6
22
6
00-03
6
4
23
6
00-05
9
2
24
6
00-10
7
4
;
run;
Data have2;
Input obs reps id $ x y;
Datalines ;
25
7
00-01
5
3
26
7
00-08
1
8
27
7
00-09
7
3
28
7
00-10
7
4
;
Run;

proc iml;
use have1 nobs nobs;
read all var {obs reps id x  y};
read all var {x y} into data;
close;
use have2;
read all var {x y} into m;
close;

center=m[:,];
start_end=t(loc(t(reps)^={.}||remove(reps,nobs)))||
          t(loc(t(reps)^=remove(reps,1)||{.}));
distance=j(nrow(reps),1,.); 

do i=1 to nrow(start_end);
 idx=start_end[i,1]:start_end[i,2];
 xx=data[idx,];
 cov=cov(xx);
 distance[idx]= mahalanobis(xx,center,cov);
end;

create want var {obs reps id x y distance};
append;
close;

quit;
SWEETSAS
Obsidian | Level 7

Thanks very much Ksharp for sharing the program. 

 

The program works, but it's computing dintance from each observation in HAVE1 to the mean of the variables in HAVE2.

 

But what I want to compute is the Mahalonobis distance from the means of the variables in each replicate in HAVE1 to the means of the variables in HAVE2, so that when the column containing distance is added back to HAVE1, it should be one value of distance for each replicate. This is because the distance is computed from the means of the variables in each replicate in HAVE1 which are unique to each replicate.

 

Thanks,

 

J

 

Ksharp
Super User
OK. No problem.
Data have1;
Input obs reps id $ x  y;
Datalines ;
1
1
00-01
5
3
2
1
00-02
4
6
3
1
00-03
6
4
4
1
00-06
4
6
5
2
00-02
4
6
6
2
00-03
6
4
7
2
00-04
7
5
8
2
00-05
9
2
9
3
00-01
5
3
10
3
00-02
4
6
11
3
00-05
9
2
12
3
00-10
7
4
13
4
00-03
6
4
14
4
00-05
9
2
15
4
00-09
7
3
16
4
00-10
7
4
17
5
00-03
6
4
18
5
00-04
7
5
19
5
00-07
5
8
20
5
00-09
7
3
21
6
00-02
4
6
22
6
00-03
6
4
23
6
00-05
9
2
24
6
00-10
7
4
;
run;
Data have2;
Input obs reps id $ x y;
Datalines ;
25
7
00-01
5
3
26
7
00-08
1
8
27
7
00-09
7
3
28
7
00-10
7
4
;
Run;

proc iml;
use have1 nobs nobs;
read all var {reps};
read all var {x y} into data;
close;
use have2;
read all var {x y} into m;
close;

center=m[:,];
level=unique(reps);
want=j(ncol(level),4,.); 

do i=1 to ncol(level);
 want[i,1]=level[i];
 idx=loc(reps=level[i]);
 want[i,2:3]=data[idx,][:,];
end;
 
 xx=want[,2:3];
 cov=cov(xx);
 want[,4]= mahalanobis(xx,center,cov);

create want from want[c={reps  x y distance}];
append from want;
close;

quit;


SWEETSAS
Obsidian | Level 7

Thanks! Thanks!! Thanks!!!. It's perfect!!!!!!

 

Please, I have a few questions on some parts of the code to enable me generalise the code to handle any number of variables: This is because a data set can have any number of variable--not just x and y; So that I will use something like &varlist.

 

want=j(ncol(level),4,.); What is the 4 doing; I assume is the number of observation in each replicate, right? Or is it the 4th column?

do i=1 to ncol(level);
want[i,1]=level[i];
idx=loc(reps=level[i]);
want[i,2:3]=data[idx,][:,]; What is the 2:3 doing? is it indexing column 2 to 3? If yes, in which data set?
end;

xx=want[,2:3]; what is the 2:3 doing?
cov=cov(xx);

 

Thanks for you incredible help,

 

J
want[,4]= mahalanobis(xx,center,cov); What is the 4 doing?

 

Thanks very much for your help.

Ksharp
Super User
OK. This code you don't need to consider about the number of variables.
Just define the calculated variables(x y ........ ) at :

read all var {x y .........} into data[c=vnames];

and Here :
read all var {x y ........} into m;





Data have1;
Input obs reps id $ x  y;
Datalines ;
1
1
00-01
5
3
2
1
00-02
4
6
3
1
00-03
6
4
4
1
00-06
4
6
5
2
00-02
4
6
6
2
00-03
6
4
7
2
00-04
7
5
8
2
00-05
9
2
9
3
00-01
5
3
10
3
00-02
4
6
11
3
00-05
9
2
12
3
00-10
7
4
13
4
00-03
6
4
14
4
00-05
9
2
15
4
00-09
7
3
16
4
00-10
7
4
17
5
00-03
6
4
18
5
00-04
7
5
19
5
00-07
5
8
20
5
00-09
7
3
21
6
00-02
4
6
22
6
00-03
6
4
23
6
00-05
9
2
24
6
00-10
7
4
;
run;
Data have2;
Input obs reps id $ x y;
Datalines ;
25
7
00-01
5
3
26
7
00-08
1
8
27
7
00-09
7
3
28
7
00-10
7
4
;
Run;

proc iml;
use have1 nobs nobs;
read all var {reps};
read all var {x y} into data[c=vnames];
close;
use have2;
read all var {x y} into m;
close;

center=m[:,];
level=unique(reps);
want=j(ncol(level),ncol(data)+2,.); 

do i=1 to ncol(level);
 want[i,1]=level[i];
 idx=loc(reps=level[i]);
 want[i,2:ncol(data)+1]=data[idx,][:,];
end;
 
 xx=want[,2:ncol(data)+1];
 cov=cov(xx);
 want[,ncol(data)+2]= mahalanobis(xx,center,cov);

names={reps}||vnames||{distance};
create want from want[c=names];
append from want;
close;

quit;

SWEETSAS
Obsidian | Level 7

The program worked as expected and computes the Mahalanobis distance when the numbers are not standardized.

 

However, when the numbers are standardized it gives an error and stops. The error is:

 

"ERROR: (execution) Matrix should be non-singular."

 

I suspect it's because one should be computing Euclidian distance when the numbers are standardized since standardized variables invlove identity matrix. Does anyone know why this is the case and how this may overcome in computing distance for standardized variables? 

 

Any help is appreciated.

 

Thank,

 

J

 

 

Ksharp
Super User
I think if they were standardized variables. You don't need calculated Mahalanobis distance,  Euclidian distance is right enough.
Mahalanobis distance is about firstly standardized  multi-dimension(multi-variables) data and calculated  Euclidian distance between them.

More info Check Rick's blog about Mahalanobis distance:
http://blogs.sas.com/content/iml/2012/02/15/what-is-mahalanobis-distance.html

And can you post full LOG information?


SWEETSAS
Obsidian | Level 7

below is the log:

 

Or could it be that it's because HAVE2  is only the mean values? The actual aim is to compute distance between of the datasets in each replicate in HAVE1 to the dataset HAVE2. I have only use the mean of HAVE2 becuase I thought that the center should be sufficient. It could be that it's because we are using only the covariance of HAVE1 instead of the pooled covariance of each of the replicate in HAVE1 and HAVE2.

 

Thanks

Below is the SAS log.

 

270 proc iml;
NOTE: IML Ready
271 use stdwant nobs nobs;
272 read all var {reps};
273 read all var {x y} into data[c=vnames];
274 close;
NOTE: Closing WORK.STDWANT
275 use stdhave;
276 read all var {x y} into m;
277 close;
NOTE: Closing WORK.STDHAVE
278
279 center=m[:,];
280 level=unique(reps);
281 want=j(ncol(level),ncol(data)+2,.);
282
283 do i=1 to ncol(level);
284 want[i,1]=level[i];
285 idx=loc(reps=level[i]);
286 want[i,2:ncol(data)+1]=data[idx,][:,];
287 end;
288
289 xx=want[,2:ncol(data)+1];
290 cov=cov(xx);
291 want[,ncol(data)+2]= mahalanobis(xx,center,cov);
NOTE: Module MAHALANOBIS loaded from the storage SASHELP.IMLMLIB.
ERROR: (execution) Matrix should be non-singular.

operation : TRISOLV at offset 8 column 18
operands : *LIT1242, U, _TEM1001

*LIT1242 1 row 1 col (numeric)

2

U 2 rows 2 cols (numeric)

0 0
0 0
_TEM1001 2 rows 30 cols (numeric)

statement : ASSIGN at offset 8 column 7
traceback : module MAHALANOBIS at offset 8 column 7

NOTE: Paused in module MAHALANOBIS.
292
293 names={reps}||vnames||{distance};
294 create want from want[c=names];
ERROR: Matrix want has not been set to a value.

statement : CREATE at line 294 column 1
295 append from want;
ERROR: No data set is currently open for output.

statement : APPEND at line 295 column 1
296 close;
297
298 quit;
NOTE: Exiting IML.
NOTE: The SAS System stopped processing this step because of errors.
NOTE: PROCEDURE IML used (Total process time):
real time 0.02 seconds
cpu time 0.04 seconds

 

Ksharp
Super User
HaHa. I know what is going on after testing the example in URL you gave. Matlab's mahalanobis distance is totally different than SAS's . and even more COV() in matlab(DF is n) is different than SAS's(DF is n-1) either . I rewrite it again by IML . Check it .


proc iml;

A={2  2,
   2  5,
   6  5,
   7  3,
   4  7,
   6  4,
   5  3,
   4  6,
   2  5,
   1  3};
                  
B={6  5,
   7  4,
   8  7,
   5  6,
   5  4};
   
   
start new_cov(x);
 Xc=x-x[:,];
 c=Xc`*Xc/nrow(x);
 return (c);
finish;   
   
start new_mahalanobis(A,B);
 xDiff=A[:,]-B[:,];  
 cA=new_cov(A);
 cB=new_cov(B);
 pC=(nrow(A)/(nrow(A)+nrow(B)))#cA +
    (nrow(B)/(nrow(A)+nrow(B)))#cB ;
 d=sqrt(xDiff*inv(pC)*xDiff`); 
 return (d);
finish;
                    
 dis=new_mahalanobis(A,B);
print dis;
quit;

Ksharp
Super User
Here is what I got according to your posted URL. Enjoy!
Data have1;
Input obs reps id $ x  y;
Datalines ;
1
1
00-01
5
3
2
1
00-02
4
6
3
1
00-03
6
4
4
1
00-06
4
6
5
2
00-02
4
6
6
2
00-03
6
4
7
2
00-04
7
5
8
2
00-05
9
2
9
3
00-01
5
3
10
3
00-02
4
6
11
3
00-05
9
2
12
3
00-10
7
4
13
4
00-03
6
4
14
4
00-05
9
2
15
4
00-09
7
3
16
4
00-10
7
4
17
5
00-03
6
4
18
5
00-04
7
5
19
5
00-07
5
8
20
5
00-09
7
3
21
6
00-02
4
6
22
6
00-03
6
4
23
6
00-05
9
2
24
6
00-10
7
4
;
run;
Data have2;
Input obs reps id $ x y;
Datalines ;
25
7
00-01
5
3
26
7
00-08
1
8
27
7
00-09
7
3
28
7
00-10
7
4
;
Run;

proc iml;
use have1 nobs nobs;
read all var {reps};
read all var {x y} into data;
close;
use have2;
read all var {x y} into B;
close;

  
start new_cov(x);
 Xc=x-x[:,];
 c=Xc`*Xc/nrow(x);
 return (c);
finish;   
   
start new_mahalanobis(A,B);
 xDiff=A[:,]-B[:,];  
 cA=new_cov(A);
 cB=new_cov(B);
 pC=(nrow(A)/(nrow(A)+nrow(B)))#cA +
    (nrow(B)/(nrow(A)+nrow(B)))#cB ;
 d=sqrt(xDiff*inv(pC)*xDiff`); 
 return (d);
finish;



start_end=t(loc(t(reps)^={.}||remove(reps,nobs)))||
          t(loc(t(reps)^=remove(reps,1)||{.}));
want=j(nrow(start_end),2,.);
 
want[,1]=reps[start_end[,1]]; 
do i=1 to nrow(start_end);
 idx=start_end[i,1]:start_end[i,2];
 A=data[idx,];
 want[i,2]=new_mahalanobis(A,B);
end;

create want from want[c={reps distance}];
append from want;
close;
quit;


SWEETSAS
Obsidian | Level 7

So means of x and y in HAVE2 is the center of Mahalanobis distance ?

Yes. HAVE2 is the center of Mahalanobis distance.
SWEETSAS
Obsidian | Level 7
http://people.revoledu.com/kardi/tutorial/Similarity/MahalanobisDistance.html

The above link contains exactly what I want to do. That is to compute distance between two matrix. In my case HAVE1 and HAVE2. The only difference is that HAVE1 have multiple data set in it, so I am computing distance between each of those replicate in HAVE1 and HAVE2. In the link above, the data were standardized and that didnt cause any problem. The site provide MATLAB code. But I want it in SAS.

Thanks,

J

Rick_SAS
SAS Super FREQ

I don't think this latest issue has anything to do with whether the variables are standardized. If you standardize a variable, it merely centers data and scales the data to unit variance.  That makes the covariance matrix equal to the correlation matrix of the unstandardized data.

 

I suspect that you have degenerate data for one of the subgroups.  For example, constant data is degenerate because it has no variance.

In general, any collinear data has a singular covariance matrix and is degenerate.

 

Print the data for the group that is giving you problems. I think you'll find that it is degenerate.

 

[BTW, the Mahalanobis distance is different from the Euclidean distance unless the covariance matrix equals the identity matrix.]

 

sas-innovate-2024.png

Join us for SAS Innovate April 16-19 at the Aria in Las Vegas. Bring the team and save big with our group pricing for a limited time only.

Pre-conference courses and tutorials are filling up fast and are always a sellout. Register today to reserve your seat.

 

Register now!

Multiple Linear Regression in SAS

Learn how to run multiple linear regression models with and without interactions, presented by SAS user Alex Chaplin.

Find more tutorials on the SAS Users YouTube channel.

From The DO Loop
Want more? Visit our blog for more articles like these.
Discussion stats
  • 23 replies
  • 3091 views
  • 9 likes
  • 4 in conversation