Threaded Mode | Linear Mode

Rodger Rosenbaum · (This post was last modified: 12-12-2021 08:35 AM by Rodger Rosenbaum.)

Note that the first thing in the little program I gave in the previous post is the LR function.

If the data are entered into the statistics registers in the usual way, the advantage of doing it this way is that you can have access to the two ways of fitting. If you execute LR, you get the ordinary least squares result; if you run the ORTH program, you get orthogonal least squares; comparison of the two is easy.

MNH · (This post was last modified: 12-18-2021 10:15 PM by MNH.)

(12-11-2021 05:38 PM)Rodger Rosenbaum Wrote: The result will be the slope and intercept for the line which minimizes the SSE with the deviation from the line in the orthogonal direction.

My results are:

120.175378956
-80799.308982

which looks strange to me, probably because I don't know a lot about statistics. My ΣDAT is

y

29945.480
30002.951
30006.678
30058.926
30114.903
30119.876
30221.977

x

921.773
922.245
982.237
921.687
923.001
983.350
924.059

In the matrix editor I entered the 'x' data in column 1 and the 'y' data in column 2. Is that correct?

The actual survey data is

PT # Northing Easting Elevation Description

248 29945.480 921.773 0.000 PROP_COR
249 30002.951 922.245 0.000 PROP_COR
250 30006.678 982.237 0.000 PROP_COR
251 30058.926 921.687 0.000 PROP_COR
252 30114.903 923.001 0.000 PROP_COR
253 30119.876 983.350 0.000 PROP_COR
254 30221.977 924.059 0.000 PROP_COR

Please see the attached file (thumbnail) to see the relationship of the points. North is straight up.

MNH · 12-19-2021, 12:10 AM

(12-11-2021 05:38 PM)Rodger Rosenbaum Wrote: << LR SWAP DROP DUP CORR DUP SIGN 4 ROLLD SQ / SWAP INV -
2 / DUP SQ 1 + √ ROT * + DUP ΣX * NEG ΣY + NΣ / >>

Bytes: 98.5
Checksum: # 55427d

Are these numbers correct?

MNH · 12-19-2021, 12:31 AM

(12-11-2021 12:13 PM)Rodger Rosenbaum Wrote: The slope and intercept for the orthogonal fit are:

slope = .010414546538
intercept = 1525424.85897

The slope needs to relate to a clockwise direction from 0 - 360 degrees.

Rodger Rosenbaum · 12-19-2021, 12:46 AM

(12-19-2021 12:10 AM)MNH Wrote:
(12-11-2021 05:38 PM)Rodger Rosenbaum Wrote: << LR SWAP DROP DUP CORR DUP SIGN 4 ROLLD SQ / SWAP INV -
2 / DUP SQ 1 + √ ROT * + DUP ΣX * NEG ΣY + NΣ / >>

Bytes: 98.5
Checksum: # 55427d

Are these numbers correct?

I am using this on a HP49g+, so the byte count may be different than what you get on your simulated HP48G.

I get:

Bytes: 90.
Checksum: # E548h

I triple checked the program and what you have is correct.

Rodger Rosenbaum · (This post was last modified: 12-19-2021 01:33 AM by Rodger Rosenbaum.)

(12-19-2021 12:31 AM)MNH Wrote:
(12-11-2021 12:13 PM)Rodger Rosenbaum Wrote: The slope and intercept for the orthogonal fit are:

slope = .010414546538
intercept = 1525424.85897

The slope needs to relate to a clockwise direction from 0 - 360 degrees.

I'm using the standard analytic geometry slope/intercept form y = mx + b, where m is the slope and b is the intercept: https://www.chilimath.com/lessons/interm...cept-form/

This is also what the LR function uses.

It shouldn't be too hard to convert to what you need.

Rodger Rosenbaum · 12-19-2021, 01:08 AM

(12-18-2021 10:13 PM)MNH Wrote:
(12-11-2021 05:38 PM)Rodger Rosenbaum Wrote: The result will be the slope and intercept for the line which minimizes the SSE with the deviation from the line in the orthogonal direction.

My results are:

120.175378956
-80799.308982

which looks strange to me, probably because I don't know a lot about statistics. My ΣDAT is

y

29945.480
30002.951
30006.678
30058.926
30114.903
30119.876
30221.977

x

921.773
922.245
982.237
921.687
923.001
983.350
924.059

In the matrix editor I entered the 'x' data in column 1 and the 'y' data in column 2. Is that correct?

The actual survey data is

PT # Northing Easting Elevation Description

248 29945.480 921.773 0.000 PROP_COR
249 30002.951 922.245 0.000 PROP_COR
250 30006.678 982.237 0.000 PROP_COR
251 30058.926 921.687 0.000 PROP_COR
252 30114.903 923.001 0.000 PROP_COR
253 30119.876 983.350 0.000 PROP_COR
254 30221.977 924.059 0.000 PROP_COR

Please see the attached file (thumbnail) to see the relationship of the points. North is straight up.

After you enter the data into ΣDAT you need to get into the statistics menu and execute the "Fit data" function (it's the 3rd function in the list). Doing this will fill in all the summary statistics such as ΣX, ΣY, ΣXY, etc. Without doing this the LR function, and the ORTH program won't work.

A check for all the data being entered is this: Recall ΣDAT to the stack; you should see a two column, 7 row matrix. Execute the calculator's ABS function. This will calculate the square root of the sum of the squares (the Frobenius norm) of all the elements of the ΣDAT matrix. You should get 79589.685862

If this is ok, try the LR function (ordinary least squares) and you should get on the stack (reverse order from what ORTH will get):

Intercept: 30074.5044063
Slope: (-7.71315429794E-3)

Now clear the stack and execute ORTH; you should get on the stack:

-1161.69645637
1121788.42831

I also checked the ORTH result with a method using the singular value decomposition and got exactly the same result.

Joe Horn · (This post was last modified: 12-19-2021 03:20 AM by Joe Horn.)

(12-19-2021 12:46 AM)Rodger Rosenbaum Wrote:
(12-19-2021 12:10 AM)MNH Wrote: Bytes: 98.5
Checksum: # 55427d

Are these numbers correct?

I am using this on a HP49g+, so the byte count may be different than what you get on your simulated HP48G.

I get:

Bytes: 90.
Checksum: # E548h

I triple checked the program and what you have is correct.

Rodger, your program's numeric literals are integer-type objects, which the 48G does not have. Replacing them with reals on the 49g+/50g changes the checksum to #27D3h. On a 48G, the checksum is #D883h. All are 90 bytes.

Rodger Rosenbaum · 12-19-2021, 03:54 AM

(12-19-2021 03:19 AM)Joe Horn Wrote:
(12-19-2021 12:46 AM)Rodger Rosenbaum Wrote: I am using this on a HP49g+, so the byte count may be different than what you get on your simulated HP48G.

I get:

Bytes: 90.
Checksum: # E548h

I triple checked the program and what you have is correct.

Rodger, your program's numeric literals are integer-type objects, which the 48G does not have. Replacing them with reals on the 49g+/50g changes the checksum to #27D3h. On a 48G, the checksum is #D883h. All are 90 bytes.

Would the checksum on an Emu48 be the same as on a real 48G?

Joe Horn · 12-19-2021, 09:04 AM

(12-19-2021 03:54 AM)Rodger Rosenbaum Wrote:
(12-19-2021 03:19 AM)Joe Horn Wrote: Rodger, your program's numeric literals are integer-type objects, which the 48G does not have. Replacing them with reals on the 49g+/50g changes the checksum to #27D3h. On a 48G, the checksum is #D883h. All are 90 bytes.

Would the checksum on an Emu48 be the same as on a real 48G?

I hope somebody who uses it can answer that for you. Unfortunately, I don't, so I can't.

MNH · 12-19-2021, 03:23 PM

I appreciate all of the help being offered here! I may be a victim of misunderstanding the mathematical lexicon involved with finding a best fit line (BFL). I recently read on a website that a BFL is an orthogonal projection using least squares. Please look at the attached thumbnail in this thread. The BFL (248 going N. to 254) on the W. side of the street is

248 29945.480 921.773
249 30002.951 922.245
251 30058.926 921.687
252 30114.903 923.001
254 30221.977 924.059

The office software I'm using creates 2 points, which define the BFL as being N 00°28'36" E (azimuth 00°28'36" or clockwise 00°28'36" from the N.). The first of these points is near 248, the second is near 254. The orthogonal offsets to the BFL are

Point: 248 Offset: 0.24655 Right
Point: 249 Offset: 0.24033 Right
Point: 251 Offset: 0.78342 Left
Point: 252 Offset: 0.06476 Right
Point: 254 Offset: 0.23177 Right

I believe I need I to create a new matrix (ΣDAT) in my Emu48. I believe I entered erroneous data after examining a recent post in this thread. If my thinking is correct, the resulting slope will indicate the direction of the BFL. I have to digest all of the information presented here. Again, thank you for the help!

Rodger Rosenbaum · 12-19-2021, 03:26 PM

(12-19-2021 03:23 PM)MNH Wrote: I appreciate all of the help being offered here! I may be a victim of misunderstanding the mathematical lexicon involved with finding a best fit line (BFL). I recently read on a website that a BFL is an orthogonal projection using least squares. Please look at the attached thumbnail in this thread. The BFL (248 going N. to 254) on the W. side of the street is

248 29945.480 921.773
249 30002.951 922.245
251 30058.926 921.687
252 30114.903 923.001
254 30221.977 924.059

The office software I'm using creates 2 points, which define the BFL as being N 00°28'36" E (azimuth 00°28'36" or clockwise 00°28'36" from the N.). The first of these points is near 248, the second is near 254. The orthogonal offsets to the BFL are

Point: 248 Offset: 0.24655 Right
Point: 249 Offset: 0.24033 Right
Point: 251 Offset: 0.78342 Left
Point: 252 Offset: 0.06476 Right
Point: 254 Offset: 0.23177 Right

I believe I need I to create a new matrix (ΣDAT) in my Emu48. I believe I entered erroneous data after examining a recent post in this thread. If my thinking is correct, the resulting slope will indicate the direction of the BFL. I have to digest all of the information presented here. Again, thank you for the help!

I don't see a thumbnail.

Rodger Rosenbaum · (This post was last modified: 12-19-2021 05:26 PM by Rodger Rosenbaum.)

Look at the two images part way down the page here:

https://stats.stackexchange.com/question...es-via-pca

The first image shows ordinary least squares which is what the HP48G LR function returns. The second image shows what the ORTH program returns.

I've got the full solution for you now.

Besides entering the ORTH program, create two new variables X and Y. Put your data in those variables as single column matrices, with X and Y being reversed from what was used previously. Designate X and Y parts of your latest data like this:

............X.............Y
248 29945.480 921.773
249 30002.951 922.245
251 30058.926 921.687
252 30114.903 923.001
254 30221.977 924.059

Type in another little program called OFIT:

<< X TRN DUP 1. CON AUGMENT TRN ORTH * Y - >>

Run it all like this:

Type your data into the ΣDAT matrix (as a two column matrix), run "Fit data..." in the statistics menu, type your X and Y data in separate single column matrices and store in the X and Y variables, then run OFIT.

Here's what I get when I do all this.

[ -.246561355 ]
[ -.240335279 ]
[ .783442323 ]
[ -.064763432 ]
[ -.231782261 ]

If you don't want to see so many digits, put "5 RND" at the end of the OFIT program.

The signs may be reversed from what you want--that's an easy fix, eh?

MNH · 12-19-2021, 08:12 PM

(12-19-2021 04:00 PM)Rodger Rosenbaum Wrote: Here's what I get when I do all this.

[ -.246561355 ]
[ -.240335279 ]
[ .783442323 ]
[ -.064763432 ]
[ -.231782261 ]

If you don't want to see so many digits, put "5 RND" at the end of the OFIT program.

That's just fabulous! Actually, even 5 digits after the decimal point is an overkill. I would never need more than 3 digits.

Quote:The signs may be reversed from what you want--that's an easy fix, eh?

The signs are reversed because a negative sign would indicate a point is to the left of a line. I'll put all this code together and post my program within a week.

MNH · 12-19-2021, 10:39 PM

(12-19-2021 04:00 PM)Rodger Rosenbaum Wrote: << X TRN DUP 1. CON AUGMENT TRN ORTH * Y - >>

I can't find the AUGMENT function in the HP 48G Series AUR Manual.

Rodger Rosenbaum · (This post was last modified: 12-20-2021 12:18 AM by Rodger Rosenbaum.)

(12-19-2021 10:39 PM)MNH Wrote:
(12-19-2021 04:00 PM)Rodger Rosenbaum Wrote: << X TRN DUP 1. CON AUGMENT TRN ORTH * Y - >>

I can't find the AUGMENT function in the HP 48G Series AUR Manual.

That's because AUGMENT is part of the HP49g+/HP50g additional functions. :-)

There's an easy way to deal with this. When you type in the X data normally it would be like this:

[ 29945.480 ]
[ 30002.951 ]
[ 30058.926 ]
[ 30114.903 ]
[ 30221.977 ]

But just type in an extra column of 1s like this:

[ 29945.480 1 ]
[ 30002.951 1 ]
[ 30058.926 1 ]
[ 30114.903 1 ]
[ 30221.977 1 ]

and store that in the X variable. Then the OFIT program becomes:

<< X ORTH * Y - >>

Don't add a column of 1s when you're typing the ΣDAT data.

Albert Chan · 12-20-2021, 01:11 AM

(12-20-2021 12:15 AM)Rodger Rosenbaum Wrote: ... Then the OFIT program becomes:

<< X ORTH * Y - >>

It may make more sense to return ORTH residuals, instead of vertical residuals.
Just multiply above by correction factor 1/√(1+m^2), where m = ORTH slope

1/√(1+m^2) = 1/√(1+tan(θ)^2) = cos(θ)

Bonus:

With ORTH residuals, ordering of X,Y does not matter.
Except for sign difference, we would get the same ORTH residuals.

Rodger Rosenbaum · 12-21-2021, 12:59 AM

When I saw the results of MNH's office software:

Point: 248 Offset: 0.24655 Right
Point: 249 Offset: 0.24033 Right
Point: 251 Offset: 0.78342 Left
Point: 252 Offset: 0.06476 Right
Point: 254 Offset: 0.23177 Right

The words "Right" and "Left" led me to think that the offsets were taken in a direction parallel to the coordinate axes. I didn't notice the sentence just above that said "The orthogonal offsets to the BFL are". :-(

But it's easy to adjust the offsets to be orthogonal offsets with a few additions to OFIT. Here's the new OFIT:

<< X ORTH DUP 1. GET UNROT * Y - SWAP ATAN COS * 5. RND >>

With this version I get:

[ -.24655 ]
[ -.24033 ]
[ .78342 ]
[ -.06476 ]
[ -.23177 ]

Which is exactly what MNH's office software gets, evidence that orthogonal
offsets are what they want. :-)

The variance of these data is very small and because of this almost any method of regression gets results very close to any of the others.

The result for the slope returned by ORTH is .008321172, but calculating it with Mathematica and 50 digit precision gets .00832117201288.

Calculating this problem with the summary statistics derived by the HP48/49g+/50g is not the best way, but even so the result from the ORTH and OFIT routines is returning 8 accurate digits as determined by the Mathematica result.

MNH · (This post was last modified: 12-21-2021 02:15 AM by MNH.)

(12-19-2021 04:00 PM)Rodger Rosenbaum Wrote: Type your data into the ΣDAT matrix (as a two column matrix) ...

My data stored in 'ΣDAT'

[[29945.480 921.773]
[30002.951 922.245]
[30058.926 921.687]
[30114.903 923.001]
[30221.977 924.059]]

Quote:... run "Fit data..." in the statistics menu ...

I get

'672.349+0.008*X'
Correlation: 0.894
Covariance: 94.090

Quote: type your X and Y data in separate single column matrices and store in the X and Y variables ...

'x' data stored in 'X'

[[29945.480 1]
[30002.951 1]
[30058.926 1]
[30114.903 1]
[30221.977 1]]

'y' data stored in 'Y'

[[921.773]
[922.245]
[921.687]
[923.001]
[924.059]]

Quote:... then run OFIT.

I get

- Error:
Bad Argument Type

Oops! I just realized that you posted something this morning. Executing the new OFIT program, I get

GET Error:
Bad Argument Type
4: 0.008
3: 672.345
2: 672.345
1: 1.000

MNH · 12-21-2021, 02:23 AM

(12-21-2021 01:56 AM)Rodger Rosenbaum Wrote: What do you get when you run ORTH by itself?

Sorry Rodger, I just deleted my Emu48 file. I don't understand why, but when I go to Open, File, the Emu48 takes me to Save (a file). Something went wrong. Let's continue tomorrow. Thank you for your help!