Is it possible to calculate R-squared on a total least squares regression?

Question

I am using the Deming function provided by Terry T. on this archived r-help thread. I am comparing two methods, so I have data that look like this:

y  x     stdy   stdx
1  1.2   0.23   0.67
2  1.8   0.05   0.89
4  7.5   1.13   0.44
... ...  ...   ...

I have done my Deming regression (also called "total least squares regression") and I get a slope and intercept. I would like to get a correlation coefficient so I've start calculating the $R^2$. I have manually entered the formula:

R2 <- function(coef,i,x,y,sdty){
    predy    <- (coef*x)+i
    stdyl    <- sum((y-predy)^2)   ### The calculated std like if it was a lm (SSres)
    Reelstdy <- sum(stdy)          ### the real stdy from the data  (SSres real)
    disty    <- sum((y-mean(y))^2) ### SS tot
    R2       <- 1-(stdyl/disty)    ### R2 formula
    R2avecstdyconnu <- 1-(Reelstdy/disty) ### R2 with the known stdy
    return(data.frame(R2, R2avecstdyconnu, stdy, Reelstdy))
}

This formula works and gives me output.

• Which of the two $R^2$s makes more sense? (I personally think of both of them as kind of biased.)
• Is there a way to get a correlation coefficient from a total least squared regression?

OUTPUT FROM THE DEMING REGRESSION:

Call:
deming(x = Data$DS, y = Data$DM, xstd = Data$SES, ystd = Data$SEM,     dfbeta = T)

               Coef  se(coef)         z            p
Intercept 0.3874572 0.2249302 3.1004680 2.806415e-10
Slope     1.2546922 0.1140142 0.8450883 4.549709e-02

   Scale= 0.7906686 
> 

Answer 1

To elaborate on whuber's answer above - Pearson will give you what you want. It determines how well y correlates with x using an approach that is independent of regression model:

$\rho_{X,Y}=\dfrac{cov(X,Y)}{\sigma_{X}\sigma_{Y}}$

gx.rma from the rgr package will do total least squares and calculate Pearson for you (or you can continue with Deming and do it manually).

require(rgr)
set.seed(3)
x<-rnorm(101,mean=2.3,sd=4.2)
x<-x+seq(0,100) 
set.seed(3)
y<-rnorm(101,mean=4.9,sd=1.9)
y<-y+seq(6,206,length=101)

rma<-gx.rma(x,y)
rma$corr
[1] 0.9922014

So, the basic answer to your question is, when doing total least squares, forget R-squared and just use Pearson. You can always square that if you want a result between 0 and 1. This will do everything you need.

Having said that, I will elaborate a little as I understand it feels like we should be able to calculate an R-squared equivalent.

First, let's try a normal sum of squares regression on the data using lm. Notice that it gives the same correlation coefficient as Pearson (after square rooting and only worrying about magnitude, obviously).

ols<-lm(y~x)
sqrt(summary(ols)$r.squared)
[1] 0.9922014

This is calculated from the lm model result using the traditional sum of squares approach

$R^{2}=1-\dfrac{S_{res}}{S_{tot}}$

So, provided you use the model given by lm, (Pearson)-squared and R-squared are equivalent.

However, if you use the model from the total sum of squares regression, and try to use the latter equation, you will get a slightly different result. That's obvious because normal and total least squares use different minimisation functions so give models with slightly different gradients and intercepts. (Remember, the first equation will still give the same result as it is looking at the data only.)

This is where I get hung up though. If the two equations give the same result when using the lm model, then surely there must be an equivalent formulation for the latter equation, but when using the total least squares model, which also gives the same result?

I had a quick play around with different approaches using the appropriate minimisation function (as has the poster here: Coefficient of determination of a orthogonal regression), but cannot find a way to do it - if there is a way.

Perhaps we are both getting hung up on the fact that Pearson and R-squared give the same result when using normal least sqaures - and there simply isn't a way to do R-squared on total least squares, which will give the same result as Pearson. But I don't know enough about this to say why not.

Answer 2

Using the package "mcr"

and using function to generate your deming regression model

yourmodel<-mcreg(x, y, ...) # you need to be familiar with the various types of deming constant SD or CV%. these can give very different results. But that's different question.

and producing a plot using the function

MCResult.plot(your model)

This displays the Pearson's production moment correlation on the plot for the model, which tells you the strength and the direction of the linear relationship between your two x,y variables, but does not give the proportion of the variation that is explained.

Hope that helps.