I found the example in NR_C302 NR_C302 is still NR2. Is there example for NR3?

BTW: In NR3, how to do multiple linear regression? Is there such example?

Thanks.

Hello,

I am not sure what you mean by "multiple linear regression". Do you mean calling the linear regression algorithm more than once in a program or do you mean multidimensional linear regression?

Here is an example program of how you could use fitab.h from NR3.


#include "nr3.h"
#include "fitab.h"

int main()
{
Doub a,val,una,unb;
Doub b,c;
Doub pause;

Doub array_x[] = { // Insert your data here };

Doub array_y[] = { // Insert your data here };

Int n = sizeof(array_x)/sizeof(array_x[0]);

VecDoub xx(n,array_x);
VecDoub yy(n,array_y);

Fitab linear(xx,yy); // Creat an instance of the struct

// Calculate the results
a=linear.a;
b=linear.b;
c=linear.chi2;
una=linear.siga;
unb=linear.sigb;

// Send the results to the screen
cout << "The value of A is " << a << endl;
cout << endl;
cout << "The value of b is " << b << endl;
cout << endl;
cout << "The equation is in the form y = A + bx " << endl;
cout << endl;
cout << "y = " << a << " + " << b << "x" << endl;
cout << endl;
cout << "The Chi Squared value is " << c << endl;
cout << endl;
cout << "Uncertanty of A " << una << endl;
cout << endl;
cout << "Uncertanty of b " << unb << endl;
cout << endl;
// Enter a value into the equation and calculate the result
cout << "Enter a value for x " << endl;
cin >> val;
cout << "Result is " << a+b*val << endl;

cin >> pause; // Dummy variable to keep the output screen visible
}

Not sure if that is what your after, but I thought I would try to help you out.

Thanks
Matt

Thank you very much for your reply. I mean multidimensional linear regression.

Y = a + b1*X1 + b2*X2 + ... + bp*Xp

I think your example is about Y = a + b1*X. I need to calculate standard error of multiple regression coefficients and t-ration, and so on.

Thanks.

BTW: I also should call the linear regression algorithm more than once in a program.

Why I cannot find the example for NR3 in my source code?

Hello,

I am not sure what you mean by "multiple linear regression". Do you mean calling the linear regression algorithm more than once in a program or do you mean multidimensional linear regression?

Here is an example program of how you could use fitab.h from NR3.


#include "nr3.h"
#include "fitab.h"

int main()
{
Doub a,val,una,unb;
Doub b,c;
Doub pause;

Doub array_x[] = { // Insert your data here };

Doub array_y[] = { // Insert your data here };

Int n = sizeof(array_x)/sizeof(array_x[0]);

VecDoub xx(n,array_x);
VecDoub yy(n,array_y);

Fitab linear(xx,yy); // Creat an instance of the struct

// Calculate the results
a=linear.a;
b=linear.b;
c=linear.chi2;
una=linear.siga;
unb=linear.sigb;

// Send the results to the screen
cout << "The value of A is " << a << endl;
cout << endl;
cout << "The value of b is " << b << endl;
cout << endl;
cout << "The equation is in the form y = A + bx " << endl;
cout << endl;
cout << "y = " << a << " + " << b << "x" << endl;
cout << endl;
cout << "The Chi Squared value is " << c << endl;
cout << endl;
cout << "Uncertanty of A " << una << endl;
cout << endl;
cout << "Uncertanty of b " << unb << endl;
cout << endl;
// Enter a value into the equation and calculate the result
cout << "Enter a value for x " << endl;
cin >> val;
cout << "Result is " << a+b*val << endl;

cin >> pause; // Dummy variable to keep the output screen visible
}

Not sure if that is what your after, but I thought I would try to help you out.

Thanks
Matt

Why I cannot find the example for NR3 in my source code?

I just made that example up. It is not in the source code.

Thanks
Matt

Thanks.

Are there examples for NR3?

Are there examples for NR3? The subject comes up from time to time, but as far as I know, no examples have been published by the authors other than snippets, here and there, in the text itself.

Here's an example of multiple regression using Fitsvd.

Comments in the code give a link to the site where I found a simple example problem.

The input data file and the program output are described in comments at the end of the following program file.


//
// Demo of use of fitsvd for multiple regression
// The problem is from an example on
// http://mtsu32.mtsu.edu:11308/regression/level3/multireg/example.htm
//
// The independent variables are the temperature and the time taken
// to mow a lawn, and the dependent variable is the amount of water
// consumed by the person doing the mowing.
//
// The data file format is the following:
//
// The first five lines are comments and are ignored by the program.
// The next line has two ints:
// The first int tells the number of lines with data values.
// The second int tells the number of variables (independent + dependent).
// The following lines have the values of the variables of the data points.
//
// The data file name can be given on the command line. If
// no file name is given, the program will use a file named
// "indata.txt" in the current directory.
//
// davekw7x March, 2010

#include "../code/nr3.h"
#include "../code/svd.h"
#include "../code/fitsvd.h"


//
// A simple linear combination of the vector variables
//
VecDoub funks(VecDoub_I & x)
{
VecDoub p(x.size() + 1);

p[0] = 1.0;
for (int i = 0; i < x.size(); i++) {
p[i+1] = x[i];
}
return p;
}


int main(int argc, char **argv)
{
string line;
int i;
const int num_comments = 5; // First five lines are thrown away
int num_variables, num_points;
char *inname = "indata.txt";
if (argc > 1) {
inname = argv[1];
}

ifstream infile(inname);
if (!infile) {
cerr << "Can't open file " << inname << " for reading." << endl;
return EXIT_FAILURE;
}
cout << "Opened " << inname << " for reading." << endl;
// The first five lines are comments
for(i = 0; i < num_comments; i++) {
getline(infile, line);
if (!infile) {
cerr << "Couldn't read comment line " << i+1
<< " from " << inname << endl;
return EXIT_FAILURE;
}
}
cout << line << endl;
infile >> num_points >> num_variables;
getline(infile, line); // Eat the rest of the line
if (!infile || (num_variables < 2) || (num_points < 2)) {
cout << "There was a problem reading the number of variables and points."
<< endl;
return EXIT_FAILURE;
}
cout << "Number of variables = " << num_variables
<< ", number of points = " << num_points
<< endl;

MatDoub x(num_points,num_variables-1);
VecDoub y(num_points);
VecDoub ssig(num_points);
for (i = 0; i < num_points; i++) {
for (int j = 0; j < num_variables-1; j++) {
infile >> x[i][j];
}
infile >> y[i];
if (!infile) {
cerr << "There was a problem reading data point " << i+1 << endl;
return EXIT_FAILURE;
}
ssig[i] = 0.1*y[i];
getline(infile, line); // Eat the rest of the line
}
cout << "Total number of data points read = " << i << endl;
infile.close();

Fitsvd fitsvd(x, y, ssig, funks); // Instantiate the object
fitsvd.fit(); // Do the work

cout << fixed << setprecision(2);
for (i = 0; i < fitsvd.a.size(); i++) {
cout << " b" << i << " = ";
cout << setw(7) << fitsvd.a[i] << " +-";
cout << setw(7) << sqrt(fitsvd.covar[i][i]) << endl;
}

return 0;
}

// Here's the file I used for testing:

/*
Multiple regression of water consumed versus temperature and mowing time.
Reference:
http://mtsu32.mtsu.edu:11308/regression/level3/multicorrel/example.htm
Fields are Temperature (F), Time (hours), Water (ounces)

7 3 # Number of data lines and number of values on each line
75 1.85 16 # Two independent variables, dependent variable
83 1.25 20
85 1.5 25
85 1.75 27
92 1.15 32
97 1.75 48
99 1.6 48
*/

// Here's the output using GNU g++ version 4.1.2 on Centos 5.4

/*
Number of variables = 3, number of points = 7
Total number of data points read = 7
b0 = -114.98 +- 17.62
b1 = 1.44 +- 0.16
b2 = 12.10 +- 4.15
*/


In my program, I assumed measurement errors were 10% of the measured value.
Here's the output:

Total number of data points read = 7
b0 = -114.98 +- 17.62
b1 = 1.44 +- 0.16
b2 = 12.10 +- 4.15


The web site showed the following output for this example:

Water = - 122 + 1.51*Temperature + 12.5*Mowing Time


Note that my output is consistent with the "exact" answer, taking into account the uncertainties due to measurement errors.

Different assumptions for the measurement errors give different outputs. (The way it works is that you can't use 0 for ssig values, but you can certainly use smaller values than what I assumed.)

For example, I changed the ssig assignments to an estimated measurement error of 0.1 and made a minor adjustment to the formatting

ssig[i] = 0.1;


Output was

b0 = -121.655 +- 0.525
b1 = 1.512 +- 0.005
b2 = 12.532 +- 0.155



Regards,

Dave

Thank you very much for your example code. It is helpful.

There are not example for NR3, it is hard for me to know how to use it. Why are there NR2 example, but not NR3? It is strange.

The subject comes up from time to time, but as far as I know, no examples have been published by the authors other than snippets, here and there, in the text itself.

Here's an example of multiple regression using Fitsvd.

Comments in the code give a link to the site where I found a simple example problem.

The input data file and the program output are described in comments at the end of the following program file.


//
// Demo of use of fitsvd for multiple regression
// The problem is from an example on
// http://mtsu32.mtsu.edu:11308/regression/level3/multireg/example.htm
//
// The independent variables are the temperature and the time taken
// to mow a lawn, and the dependent variable is the amount of water
// consumed by the person doing the mowing.
//
// The data file format is the following:
//
// The first five lines are comments and are ignored by the program.
// The next line has two ints:
// The first int tells the number of lines with data values.
// The second int tells the number of variables (independent + dependent).
// The following lines have the values of the variables of the data points.
//
// The data file name can be given on the command line. If
// no file name is given, the program will use a file named
// "indata.txt" in the current directory.
//
// davekw7x March, 2010

#include "../code/nr3.h"
#include "../code/svd.h"
#include "../code/fitsvd.h"


//
// A simple linear combination of the vector variables
//
VecDoub funks(VecDoub_I & x)
{
VecDoub p(x.size() + 1);

p[0] = 1.0;
for (int i = 0; i < x.size(); i++) {
p[i+1] = x[i];
}
return p;
}


int main(int argc, char **argv)
{
string line;
int i;
const int num_comments = 5; // First five lines are thrown away
int num_variables, num_points;
char *inname = "indata.txt";
if (argc > 1) {
inname = argv[1];
}

ifstream infile(inname);
if (!infile) {
cerr << "Can't open file " << inname << " for reading." << endl;
return EXIT_FAILURE;
}
cout << "Opened " << inname << " for reading." << endl;
// The first five lines are comments
for(i = 0; i < num_comments; i++) {
getline(infile, line);
if (!infile) {
cerr << "Couldn't read comment line " << i+1
<< " from " << inname << endl;
return EXIT_FAILURE;
}
}
cout << line << endl;
infile >> num_points >> num_variables;
getline(infile, line); // Eat the rest of the line
if (!infile || (num_variables < 2) || (num_points < 2)) {
cout << "There was a problem reading the number of variables and points."
<< endl;
return EXIT_FAILURE;
}
cout << "Number of variables = " << num_variables
<< ", number of points = " << num_points
<< endl;

MatDoub x(num_points,num_variables-1);
VecDoub y(num_points);
VecDoub ssig(num_points);
for (i = 0; i < num_points; i++) {
for (int j = 0; j < num_variables-1; j++) {
infile >> x[i][j];
}
infile >> y[i];
if (!infile) {
cerr << "There was a problem reading data point " << i+1 << endl;
return EXIT_FAILURE;
}
ssig[i] = 0.1*y[i];
getline(infile, line); // Eat the rest of the line
}
cout << "Total number of data points read = " << i << endl;
infile.close();

Fitsvd fitsvd(x, y, ssig, funks); // Instantiate the object
fitsvd.fit(); // Do the work

cout << fixed << setprecision(2);
for (i = 0; i < fitsvd.a.size(); i++) {
cout << " b" << i << " = ";
cout << setw(7) << fitsvd.a[i] << " +-";
cout << setw(7) << sqrt(fitsvd.covar[i][i]) << endl;
}

return 0;
}

// Here's the file I used for testing:

/*
Multiple regression of water consumed versus temperature and mowing time.
Reference:
http://mtsu32.mtsu.edu:11308/regression/level3/multicorrel/example.htm
Fields are Temperature (F), Time (hours), Water (ounces)

7 3 # Number of data lines and number of values on each line
75 1.85 16 # Two independent variables, dependent variable
83 1.25 20
85 1.5 25
85 1.75 27
92 1.15 32
97 1.75 48
99 1.6 48
*/

// Here's the output using GNU g++ version 4.1.2 on Centos 5.4

/*
Number of variables = 3, number of points = 7
Total number of data points read = 7
b0 = -114.98 +- 17.62
b1 = 1.44 +- 0.16
b2 = 12.10 +- 4.15
*/


In my program, I assumed measurement errors were 10% of the measured value.
Here's the output:

Total number of data points read = 7
b0 = -114.98 +- 17.62
b1 = 1.44 +- 0.16
b2 = 12.10 +- 4.15


The web site showed the following output for this example:

Water = - 122 + 1.51*Temperature + 12.5*Mowing Time


Note that my output is consistent with the "exact" answer, taking into account the uncertainties due to measurement errors.

Different assumptions for the measurement errors give different outputs. (The way it works is that you can't use 0 for ssig values, but you can certainly use smaller values than what I assumed.)

For example, I changed the ssig assignments to an estimated measurement error of one minute and adjusted output formatting slightly

ssig[i] = 1.0 / 60.0;


Output was

b0 = -121.655 +- 0.088
b1 = 1.512 +- 0.001
b2 = 12.532 +- 0.026



Regards,

Dave

...Why are there NR2 example, but not NR3? It is strange...

Well, I'm not sure that strange is the word I would have used. In my opinion, strangeness, like evil (or like beauty, or like pornography) is in the eye of the beholder. See Footnotes.

Also, and alas, see author Bill Press's response in this thread: http://www.nr.com/forum/showthread.php?t=1875

Regards,

Dave


Footnotes:
[1]
"In the end, we self-perceiving, self-inventing, locked-in mirages are little miracles of self-reference."
---Douglas Hofstadter
I Am a Strange Loop

[2]
"Just remember, as long as I am here, you are not the strangest person in the room."
---Ally McBeal

For your example, I use R. I got this.

In your example, is it possible to get "Std. Error" column as: 6.54035, 0.06077,1.93302?
you just have sqrt(covar(i,i)).

Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) -121.65500 6.54035 -18.601 4.92e-05 ***
V1 1.51236 0.06077 24.886 1.55e-05 ***
V2 12.53168 1.93302 6.483 0.00292 **
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Residual standard error: 1.245 on 4 degrees of freedom
Multiple R-squared: 0.9937, Adjusted R-squared: 0.9905
F-statistic: 313.2 on 2 and 4 DF, p-value: 4.027e-05


The subject comes up from time to time, but as far as I know, no examples have been published by the authors other than snippets, here and there, in the text itself.

Here's an example of multiple regression using Fitsvd.

Comments in the code give a link to the site where I found a simple example problem.

The input data file and the program output are described in comments at the end of the following program file.


//
// Demo of use of fitsvd for multiple regression
// The problem is from an example on
// http://mtsu32.mtsu.edu:11308/regression/level3/multireg/example.htm
//
// The independent variables are the temperature and the time taken
// to mow a lawn, and the dependent variable is the amount of water
// consumed by the person doing the mowing.
//
// The data file format is the following:
//
// The first five lines are comments and are ignored by the program.
// The next line has two ints:
// The first int tells the number of lines with data values.
// The second int tells the number of variables (independent + dependent).
// The following lines have the values of the variables of the data points.
//
// The data file name can be given on the command line. If
// no file name is given, the program will use a file named
// "indata.txt" in the current directory.
//
// davekw7x March, 2010

#include "../code/nr3.h"
#include "../code/svd.h"
#include "../code/fitsvd.h"


//
// A simple linear combination of the vector variables
//
VecDoub funks(VecDoub_I & x)
{
VecDoub p(x.size() + 1);

p[0] = 1.0;
for (int i = 0; i < x.size(); i++) {
p[i+1] = x[i];
}
return p;
}


int main(int argc, char **argv)
{
string line;
int i;
const int num_comments = 5; // First five lines are thrown away
int num_variables, num_points;
char *inname = "indata.txt";
if (argc > 1) {
inname = argv[1];
}

ifstream infile(inname);
if (!infile) {
cerr << "Can't open file " << inname << " for reading." << endl;
return EXIT_FAILURE;
}
cout << "Opened " << inname << " for reading." << endl;
// The first five lines are comments
for(i = 0; i < num_comments; i++) {
getline(infile, line);
if (!infile) {
cerr << "Couldn't read comment line " << i+1
<< " from " << inname << endl;
return EXIT_FAILURE;
}
}
cout << line << endl;
infile >> num_points >> num_variables;
getline(infile, line); // Eat the rest of the line
if (!infile || (num_variables < 2) || (num_points < 2)) {
cout << "There was a problem reading the number of variables and points."
<< endl;
return EXIT_FAILURE;
}
cout << "Number of variables = " << num_variables
<< ", number of points = " << num_points
<< endl;

MatDoub x(num_points,num_variables-1);
VecDoub y(num_points);
VecDoub ssig(num_points);
for (i = 0; i < num_points; i++) {
for (int j = 0; j < num_variables-1; j++) {
infile >> x[i][j];
}
infile >> y[i];
if (!infile) {
cerr << "There was a problem reading data point " << i+1 << endl;
return EXIT_FAILURE;
}
ssig[i] = 0.1*y[i];
getline(infile, line); // Eat the rest of the line
}
cout << "Total number of data points read = " << i << endl;
infile.close();

Fitsvd fitsvd(x, y, ssig, funks); // Instantiate the object
fitsvd.fit(); // Do the work

cout << fixed << setprecision(2);
for (i = 0; i < fitsvd.a.size(); i++) {
cout << " b" << i << " = ";
cout << setw(7) << fitsvd.a[i] << " +-";
cout << setw(7) << sqrt(fitsvd.covar[i][i]) << endl;
}

return 0;
}

// Here's the file I used for testing:

/*
Multiple regression of water consumed versus temperature and mowing time.
Reference:
http://mtsu32.mtsu.edu:11308/regression/level3/multicorrel/example.htm
Fields are Temperature (F), Time (hours), Water (ounces)

7 3 # Number of data lines and number of values on each line
75 1.85 16 # Two independent variables, dependent variable
83 1.25 20
85 1.5 25
85 1.75 27
92 1.15 32
97 1.75 48
99 1.6 48
*/

// Here's the output using GNU g++ version 4.1.2 on Centos 5.4

/*
Number of variables = 3, number of points = 7
Total number of data points read = 7
b0 = -114.98 +- 17.62
b1 = 1.44 +- 0.16
b2 = 12.10 +- 4.15
*/


In my program, I assumed measurement errors were 10% of the measured value.
Here's the output:

Total number of data points read = 7
b0 = -114.98 +- 17.62
b1 = 1.44 +- 0.16
b2 = 12.10 +- 4.15


The web site showed the following output for this example:

Water = - 122 + 1.51*Temperature + 12.5*Mowing Time


Note that my output is consistent with the "exact" answer, taking into account the uncertainties due to measurement errors.

Different assumptions for the measurement errors give different outputs. (The way it works is that you can't use 0 for ssig values, but you can certainly use smaller values than what I assumed.)

For example, I changed the ssig assignments to an estimated measurement error of 0.1 and made a minor adjustment to the formatting

ssig[i] = 0.1;


Output was

b0 = -121.655 +- 0.525
b1 = 1.512 +- 0.005
b2 = 12.532 +- 0.155



Regards,

Dave

...
you just have sqrt(covar(i,i))....My example was offered in response to your initial question about linear regression using NR3 functions. The Fitsvd methodology gives us the coefficients obtained by minimizing chi-squared, and it gives us the covariance matrix and it gives us the value of chi-squared.

The value of covar[i][i] is equal to the variance in the estimate of parameter b[i].

Regards,

Dave

It can out put the standard error of b[i] and F-value?

My example was offered in response to your initial question about linear regression using NR3 functions. The Fitsvd methodology gives us the coefficients obtained by minimizing chi-squared, and it gives us the covariance matrix and it goves is the value of chi-squared.

The value of covar[i][i] is equal to the variance in the estimate of parameter b[i].

Regards,

Dave

It can out put the standard error of b[i] and F-value?Look at the code for Fitsvd.
I hate to repeat myself but

The Fitsvd methodology gives us the coefficients obtained by minimizing chi-squared, and it gives us the covariance matrix and it gives us the value of chi-squared.


Period.


Regards,

Dave