Ok, the code is running.
There is something I don't understand.
Here is some input data (again x,y,weight)
27,1.124793199, 0.05076084028870085627,1.124793199, 0.050760840288700856
28,1.109467023, 0.2567362070012724
29,1.125426721, 0.5381833400915066
30,1.120063104, 0.7819498266908333
31,1.124317802, 0.9316484397540968
32,1.125225209, 0.9912791191131936
33,1.139629547, 1.0
34,1.146391211, 0.9912791191131936
35,1.152133707, 0.9316484397540968
36,1.15958384 , 0.7819498266908333
37,1.145192155, 0.5381833400915066
38,1.15645546 , 0.2567362070012724
39,1.165596309, 0.050760840288700856
When I run the code on this and call fitter.theoreticalValue(33),
it returns 1.139629547, the value that was input.
Is this correct?
Luc Maisonobe wrote:
[EMAIL PROTECTED] wrote:
Thanks for the reply and my apologies for omitting the [math] marker.
Afa the model goes, I'm not sure how to answer. What I am doing is
smoothing a curve using the loess function, and the last step is to use a
weighted least square regression on each point and its neighborhood.
OK. I consider then that you already have selected the neighborhood of
your point and computed the local weights for i and want to compute only
the quadratically fitted value for the current point. Here is a simple
way to do that (I have only added more print statements, and two
different samples, one perfect and one with your data).
The trick (which I should put in the documentation) is that often the
simplest way is to have an implementation of EstimationProblem which has
all its parameters as attributes and that build measurements using an
internal class during problem setup. This allow each measurement to
retrieve the values of the parameters and apply the model as needed to
compute theoretical values which will be compared by the solver with the
measured value they already hold.
Hope this helps.
Luc
import org.apache.commons.math.estimation.EstimatedParameter;
import org.apache.commons.math.estimation.EstimationException;
import org.apache.commons.math.estimation.LevenbergMarquardtEstimator;
import org.apache.commons.math.estimation.SimpleEstimationProblem;
import org.apache.commons.math.estimation.WeightedMeasurement;
public class QuadraticFitter extends SimpleEstimationProblem {
public static void main(String[] args) {
try {
// create the uninitialized fitting problem
QuadraticFitter fitter = new QuadraticFitter();
// // register the sample points (perfect quadratic case)
// fitter.addPoint(1.0, 5.0, 0.2);
// fitter.addPoint(2.0, 10.0, 0.4);
// fitter.addPoint(3.0, 17.0, 1.0);
// fitter.addPoint(4.0, 26.0, 0.8);
// fitter.addPoint(5.0, 37.0, 0.3);
// register the sample points (practical case)
fitter.addPoint(1.0, 1.0, 0.2);
fitter.addPoint(2.0, 3.0, 0.4);
fitter.addPoint(3.0, 2.0, 1.0);
fitter.addPoint(4.0, 6.0, 0.8);
fitter.addPoint(5.0, 4.0, 0.3);
double centerX = 3.0;
// solve the problem, using a Levenberg-Marquardt algorithm
// with default settings
LevenbergMarquardtEstimator estimator =
new LevenbergMarquardtEstimator();
estimator.estimate(fitter);
System.out.println("RMS after fitting = "
+ estimator.getRMS(fitter));
// display the results
System.out.println("a = " + fitter.getA()
+ ", b = " + fitter.getB()
+ ", c = " + fitter.getC());
System.out.println("y (" + centerX + ") = "
+ fitter.theoreticalValue(centerX));
} catch (EstimationException ee) {
System.err.println(ee.getMessage());
}
}
/** Fitter for a quadratic model to a sample of 2D points.
* <p>The model is y(x) = a x<sup>2</sup> + b x + c
* its three parameters of the model are a, b and c.</p>
*/
public QuadraticFitter() {
// three parameters of the model
a = new EstimatedParameter("a", 0.0);
b = new EstimatedParameter("b", 0.0);
c = new EstimatedParameter("c", 0.0);
// provide the parameters to the base class which
// implements the getAllParameters and
// getUnboundParameters methods
addParameter(a);
addParameter(b);
addParameter(c);
}
/** Add a sample point.
* @param x abscissa
* @param y ordinate
* @param w weight
*/
public void addPoint(double x, double y, double w) {
addMeasurement(new LocalMeasurement(x, y, w));
}
/** Get the value of the quadratic coefficient.
* @return the value of a for the quadratic model
* y = a x<sup>2</sup> + b x + c
*/
public double getA() {
return a.getEstimate();
}
/** Get the value of the linear coefficient.
* @return the value of b for the quadratic model
* y = a x<sup>2</sup> + b x + c
*/
public double getB() {
return b.getEstimate();
}
/** Get the value of the constant coefficient.
* @return the value of ac for the quadratic model
* y = a x<sup>2</sup> + b x + c
*/
public double getC() {
return c.getEstimate();
}
/** Get the theoretical value of the model for some x.
* <p>The theoretical value is the value computed using
* the current state of the problem parameters.</p>
* @param x explanatory variable
* @return the theoretical value y = a x<sup>2</sup> + b x + c
*/
public double theoreticalValue(double x) {
return ((a.getEstimate() * x + b.getEstimate()) * x
+ c.getEstimate());
}
/** Get the partial derivative of the theoretical value
* of the model for some x.
* <p>The derivative is computed using
* the current state of the problem parameters.</p>
* @param x explanatory variable
* @param parameter estimated parameter (either a, b, or c)
* @return the partial derivative dy/dp
*/
private double partial(double x, EstimatedParameter parameter) {
// since we know the only parameters are a, b and c in this
// class we simply use "==" for efficiency
if (parameter == a) {
return x * x;
} else if (parameter == b) {
return x;
} else {
return 1.0;
}
}
/** Internal measurements class.
* <p>The measurement is the y value for a fixed specified x.</p>
*/
private class LocalMeasurement extends WeightedMeasurement {
public LocalMeasurement(double x, double y, double w) {
super(w, y);
this.x = x;
}
public double getTheoreticalValue() {
// the value is provided by the model for the local x
return theoreticalValue(x);
}
public double getPartial(EstimatedParameter parameter) {
// the value is provided by the model for the local x
return partial(x, parameter);
}
private final double x;
}
private EstimatedParameter a;
private EstimatedParameter b;
private EstimatedParameter c;
}
Hi,
First of all, I have added a [math] marker on the subject line. This
list
is shared among all commons projects and this type of markers help
people
filter the messages.
I will send a usage example on the list in a few hours (late evening,
european time), when I'm back home. Would you like to have anything
special in this example ? For example what kind of model do you want to
be fitted to the x,ydata ?
Luc
Selon [EMAIL PROTECTED]:
Can anybody show me an example of a weighted least squares regression
using classes like EstimationProblem, WeightedMeasurement from
apache.commons.math?
I have data that looks like this: (x,y,weight), e.g.
1,1,0.2
2,3, 0.4
3,2, 1.0
4,6, 0.8
5,4, 0.3
---------------------------------------------------------------------
To unsubscribe, e-mail: [EMAIL PROTECTED]
For additional commands, e-mail: [EMAIL PROTECTED]
---------------------------------------------------------------------
To unsubscribe, e-mail: [EMAIL PROTECTED]
For additional commands, e-mail: [EMAIL PROTECTED]
---------------------------------------------------------------------
To unsubscribe, e-mail: [EMAIL PROTECTED]
For additional commands, e-mail: [EMAIL PROTECTED]
