Estimating Coefficients for Piecewise Function I'm interested in fitting a piecewise defined function to the following data, so that I get something like:

x = c(2013,2012,2011,2010,2009,2008,2007,2006,2005,2003,2002,2001,2000,1999,1998,1996)
y = c(397,332,377,425,467,611,882,906,1032,955,843,895,873,943,898,476)

$f(x) = \begin{cases}
d & if \mbox{ x} < 2007 \\
ax^2 + bx + c & if \mbox{ x} \geq 2007
\end{cases}
$
I also want $f(x)$ to be continuous and differentiable, which it will be except for at $x=2007$, unless we apply some constraints.
To be continuous
$
\displaystyle \lim_{x \to 2007^+} f(x) = \lim_{x \to 2007^-}= f(2007) \\
$
In my example that means
$
d = a(2007)^2+b(2007)+c
$
To be differentiable, the derivatives must equal at the knot $x=2007$.
$
0 = 2a(2007) + b \\
b = -4014a
$
I've tried to program this in R
f = function(theta,y,x){
  a = theta[1] 
  b = -4014*a
  c = theta[2]
  d = a*2007^2+b*2007+c
  mu = ifelse(x>=2007,a*x^2+b*x+c,d)
  nll = -sum(dnorm(y,mu,1,log=TRUE))
  nll
}

fit = optim(par=c(4.032e-02,1.632e+05),f,y=y,x=x,method='BFGS')

I may have bad starting values, as the objective functions's values doesn't change.  Am I going about this the right way?  Perhaps there's something more stable?
 A: I would use OLS for this problem as follows
OLS regression
m<-lm(y~I(1*(x<2007))+I(x^2*(x>=2007))+I(x*(x>=2007))+I(1*(x>=2007))+0)

In order to satisfy your constraints on differentiability, however, you can use nls
m2<-nls(y~(a*2007^2-4014*a*2007+c)*I((1)*(x<2007))+
  a*I((x^2)*(x>=2007))+(-4014*a)*I((x)*(x>=2007))+c*I((1)*(x>=2007)),
  start=list(a=2.575e+01,c=1.042e+08),control=list(tol=1e-02))

Given the constraints and the data, however, the fit is poor
plot(y~x)
lines(x,fitted.values(m),col="red")
lines(x,fitted.values(m2),col="blue")


A: Here's an answer to your question. But the answer comes from basic concepts of splines. I think looking further into that topic will be very helpful for you!
If you want to make sure that your estimated $f(x)$ is continuous and follows the form that you propose, the easiest way to do this is to expand the $x$ basis. For example, if you want $f(x)$ to be a constant if $x < 2007$ but quadratic after that, you can simply have the covariates
$x_{0,i} = 1$
$x_{1,i} = x - 2007$ if $x >= 2007$ else 0
$x_{2,i} = (x - 2007)^2$ if $x >= 2007$ else 0
Then you can use simple linear regression to fit $y \sim~ x_0 + x_1 + x_2$. To get your given predictions for a new value of $x$, you just need to get your expansion into $x_0$, $x_1$ and $x_2$ and plug your regression estimates. Note that this will insure a fit that is always continuous, constant before 2007 and quadratic after 2007. 
While I do believe this answers your question as asked, I will also comment that this is not how I would attempt to model a non-linear effect; you've built very strong restrictions into your model, which I'm guessing were motivated post-hoc by examining your data. This would make the generalization of your results very questionable.  Using standard spline methods would lead to more reliable estimation. 
