Sales<-read.csv("iPad.csv")

Sales
##     Sales
## 1   0.270
## 2   1.119
## 3   2.315
## 4   1.783
## 5   3.717
## 6   6.892
## 7   4.373
## 8   3.793
## 9   5.288
## 10  7.347
## 11  8.737
## 12  9.752
## 13  9.998
## 14 14.142
## 15 15.535
## 16 19.647
## 17 20.338
## 18 17.073
## 19 27.044
## 20 35.064
## 21 24.428
## 22 26.919
## 23 47.189
Sales=ts(Sales,start=c(2014,3),freq=4)
#Plot data
library(ggplot2) 
plot(Sales,type="l",lty=2, col="red", ylab="", xlab="")
points(Sales,pch=20, col="blue")
title("Quarterly sales (millions)")

#Plot cumulative sales

Y=cumsum(Sales)
Y=ts(Y,start=c(2014,3),freq=4)
plot(Y,type="l", lty=2, col="red", ylab="",xlab="")
points(Y,pch=20,col="blue")
title("Cumulative sales (millions)")

#Bass model and compute m,p,q
Y=c(0,Y[1:(length(Y)-1)])
Ysq=Y**2
out=lm(Sales~Y+Ysq)
summary(out)
## 
## Call:
## lm(formula = Sales ~ Y + Ysq)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -8.4210 -1.0372  0.3226  1.3471  8.8555 
## 
## Coefficients:
##               Estimate Std. Error t value Pr(>|t|)    
## (Intercept)  2.1630612  1.4464707   1.495 0.150421    
## Y            0.1640106  0.0363801   4.508 0.000215 ***
## Ysq         -0.0001047  0.0001488  -0.704 0.489611    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 4.062 on 20 degrees of freedom
## Multiple R-squared:  0.8984, Adjusted R-squared:  0.8882 
## F-statistic:  88.4 on 2 and 20 DF,  p-value: 1.175e-10
a=out$coef[1]
b=out$coef[2]
c=out$coef[3]
mplus=(-b+sqrt(b**2-4*a*c))/(2*c)
mminus=(-b-sqrt(b**2-4*a*c))/(2*c)
m=mminus
p=1/m
q=b+p
#Compute Bass Model
bassModel=function(p,q,m,T=100)
{
  S=double(T)
  Y=double(T+1)
  Y[1]=0
  for(t in 1:T)
  {
    S[t]=p*m+(q-p)*Y[t]-(q/m)*Y[t]**2
  Y[t+1]=Y[t]+S[t]
  }
  return(list(sales=S,cumSales=cumsum(S)))
}

Spred=bassModel(p,q,m,T=23)$sales

Spred=ts(Spred,start=c(2014,3),freq=4)
ts.plot(Sales,Spred,col=c("blue","red"))
legend("topleft", legend=c("Actual", "Bass Model"),fill=c("blue","red"))

#Proyection

Spred=bassModel(p,q,m)$sales
CumSpred=ts(cumsum(Spred),start=c(2014,3),freq=4)
CumSales=ts(cumsum(Sales),start=c(2014,3),freq=4)
ts.plot(CumSales,CumSpred,col=c("blue","red"))
legend("topleft", legend=c("Actual", "Bass Model"),fill=c("blue","red"))
title("Predicted Cumulative Sales IPads")