Note that 
is more often written as 
or just 
.
In the case where 
is unknown and you want to estimate it, you can use nonlinear least squares on this.
[In Matlab, see lsqcurvefit, for example]
Alternatively, since given a value for 
you can write it as a linear regression problem, you can write the problem as partially linear model*, where for example, given some value of 
you can estimate 
and 
by least squares, so you can just optimize the MSE over 
). Any univariate optimizer should work nicely for that. While 
and 
can actually be eliminated, you can avoid the effort of doing that: within the function that calculates the sum of squares of residuals for a given 
, you just compute the least squares fit of 
and 
in a regression of 
on 
and the sum of squares of that fit is the function value.
[In Matlab, see fminsearch as an example of an optimizer, but since it's sum of squares, you should be able to take advantage of lsqnonlin]
* not to be confused with partial least squares which is quite a different thing.
I don't have Matlab at present, but I can illustrate some of these ideas in R easily enough; it's not hard to do the same kind of thing in Matlab.
First I generated some data with 
,
and 
:
set.seed(329783)
m=11; a=0.75; b=5
x = runif(100,3,21)
y = a*abs(x-m)+b + rnorm(100,0,.3)
a) If you have moderately good guesses at the parameters, this is simple with a nonlinear least squares program:
> Vfit0 = nls(y~a*abs(x-m)+b,start=list(m=10,a=1,b=4))
> summary(Vfit0)
Formula: y ~ a * abs(x - m) + b
Parameters:
Estimate Std. Error t value Pr(>|t|)
m 11.03212 0.03891 283.5 <2e-16 ***
a 0.74747 0.01122 66.6 <2e-16 ***
b 5.01031 0.05792 86.5 <2e-16 ***
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
Residual standard error: 0.2827 on 97 degrees of freedom
Number of iterations to convergence: 4
Achieved convergence tolerance: 2.878e-09
If you don't have good guesses at the parameters (or you want to completely automate it), good estimates are easy enough to construct for this problem.
If you subtract off the minimum 
(to almost get rid of 
), the square of 
is approximately quadratic in 
(it's not all that sensitive to accuracy in 
). The x-value for the minimum of that quadratic should still be near 
.
The minimum of a quadratic ($a_2x^2+a_1x+a_0$) will be at 
.
So calculating that for a quadratic fit (in 
) to 
, should give a good estimate of 
. Then given that 
, a linear regression of 
on 
gives a good 
and 
. So those should be nice starting values --
Here's how I did that in R:
qmin = function(v) -v[2]/v[3]/2
y0 = y^2-min(y^2)
m0 = qmin( lm( y0~x+I(x^2) )$coefficients ) # start estimate of m
stcoef = lm(y~abs(x-m0))$coefficients # to get start estimates of a and b
Vfit = nls(y~sign(x-m)*a*(x-m)+b,start=list(m=m0,a=stcoef[2],b=stcoef[1]))
the final fit was the same as above so I won't give that, but the starting values were
m a b
11.14097 0.7536526 4.984033
which are much closer to the final estimates than my (quite sufficient) original "guess"
Answer from Glen_b on Stack ExchangeVideos
Note that is more often written as or just .
In the case where is unknown and you want to estimate it, you can use nonlinear least squares on this.
[In Matlab, see lsqcurvefit, for example]
Alternatively, since given a value for you can write it as a linear regression problem, you can write the problem as partially linear model*, where for example, given some value of you can estimate and by least squares, so you can just optimize the MSE over ). Any univariate optimizer should work nicely for that. While and can actually be eliminated, you can avoid the effort of doing that: within the function that calculates the sum of squares of residuals for a given , you just compute the least squares fit of and in a regression of on and the sum of squares of that fit is the function value.
[In Matlab, see fminsearch as an example of an optimizer, but since it's sum of squares, you should be able to take advantage of lsqnonlin]
* not to be confused with partial least squares which is quite a different thing.
I don't have Matlab at present, but I can illustrate some of these ideas in R easily enough; it's not hard to do the same kind of thing in Matlab.
First I generated some data with , and :
set.seed(329783)
m=11; a=0.75; b=5
x = runif(100,3,21)
y = a*abs(x-m)+b + rnorm(100,0,.3)
a) If you have moderately good guesses at the parameters, this is simple with a nonlinear least squares program:
> Vfit0 = nls(y~a*abs(x-m)+b,start=list(m=10,a=1,b=4))
> summary(Vfit0)
Formula: y ~ a * abs(x - m) + b
Parameters:
Estimate Std. Error t value Pr(>|t|)
m 11.03212 0.03891 283.5 <2e-16 ***
a 0.74747 0.01122 66.6 <2e-16 ***
b 5.01031 0.05792 86.5 <2e-16 ***
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
Residual standard error: 0.2827 on 97 degrees of freedom
Number of iterations to convergence: 4
Achieved convergence tolerance: 2.878e-09
If you don't have good guesses at the parameters (or you want to completely automate it), good estimates are easy enough to construct for this problem.
If you subtract off the minimum (to almost get rid of ), the square of is approximately quadratic in (it's not all that sensitive to accuracy in ). The x-value for the minimum of that quadratic should still be near .
The minimum of a quadratic ($a_2x^2+a_1x+a_0$) will be at .
So calculating that for a quadratic fit (in ) to , should give a good estimate of . Then given that , a linear regression of on gives a good and . So those should be nice starting values --
Here's how I did that in R:
qmin = function(v) -v[2]/v[3]/2
y0 = y^2-min(y^2)
m0 = qmin( lm( y0~x+I(x^2) )$coefficients ) # start estimate of m
stcoef = lm(y~abs(x-m0))$coefficients # to get start estimates of a and b
Vfit = nls(y~sign(x-m)*a*(x-m)+b,start=list(m=m0,a=stcoef[2],b=stcoef[1]))
the final fit was the same as above so I won't give that, but the starting values were
m a b
11.14097 0.7536526 4.984033
which are much closer to the final estimates than my (quite sufficient) original "guess"
Answer from Glen_b on Stack ExchangeAs mentioned by others, the least-squares problem is much easier to solve. But there’s another important reason: assuming IID Gaussian noise, the least-squares solution is the Maximum-Likelihood estimate.
$$\min_{a,b}\sum_{k=1}^n(ax_k+b-y_k)^2$$ has a simple analytical solution.
$$\min_{a,b}\sum_{k=1}^n|ax_k+b-y_k|$$ is difficult.
One of reasons is that the absolute value is not differentiable.