Generación de datos

# Generación de la curva de demanda (linealmente dependiente del precio)
demanda = function(p, alpha = -40, beta = 500, sd = 10) {
  error = rnorm(length(p), sd = sd)
  q = p*alpha + beta + error
}

Generación de los datos a analizar

set.seed(100)
 hist.precio = rnorm(252, mean = 6, sd = .5)
 hist.demanda = demanda(hist.precio)
 hist.ingresos = hist.precio*hist.demanda
 costes_unitarios = 4
 hist.costes = costes_unitarios*hist.demanda
 hist.beneficios = (hist.precio - costes_unitarios)*hist.demanda # From the precio equation

df = data.frame('periodo' = seq(1,252),'diario.precio' = hist.precio,
 'demanda_diaria' = hist.demanda, 'ingresos_diarios' = hist.ingresos,
 'costes_diarios' = hist.costes, 'beneficios_diarios' = hist.beneficios)

library(ggplot2)
 ggplot(df, aes(periodo, ingresos_diarios, colour = 'ingresos')) + theme_minimal() + geom_line() + geom_line(aes(periodo, beneficios_diarios, colour = 'beneficios')) + geom_line(aes(periodo, costes_diarios, colour = 'costes')) + labs(title = 'Historico de costes, beneficios e ingresos', colour = '')

Generación del modelo

modelo_lineal = lm(hist.demanda ~ hist.precio) # linear model for demanda
summary(modelo_lineal)

## 
## Call:
## lm(formula = hist.demanda ~ hist.precio)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -32.990  -6.100   0.503   7.472  33.837 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)    
## (Intercept)  489.177      8.287   59.03   <2e-16 ***
## hist.precio  -38.314      1.379  -27.79   <2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 10.5 on 250 degrees of freedom
## Multiple R-squared:  0.7555, Adjusted R-squared:  0.7545 
## F-statistic: 772.4 on 1 and 250 DF,  p-value: < 2.2e-16

# Estimación de parametros
beta = modelo_lineal$coefficients[1]
alpha = modelo_lineal$coefficients[2]  

# Maximización del beneficio (Derivada del  beneficio con respecto al precio)
p.ingresos = -beta/(2*alpha) 

# Maximización de los ingresos (Derivada de los ingresos con respecto al precio)
p.beneficios = (alpha*costes_unitarios - beta)/(2*alpha) # Estimación de beneficios

# Representación de como evolucionan los ingresos con respecto al precio
modelado_ingresos = function(p) p*(modelo_lineal$coefficients[2]*p + modelo_lineal$coefficients[1])

# Representación de como evolucionan los ingresos con respecto al beneficio
modelado_beneficios = function(p) (p - costes_unitarios)*(modelo_lineal$coefficients[2]*p + modelo_lineal$coefficients[1])

# Estimación del punto donde se maximizan los ingresos para el precio óptimo calculado
max_ingresos = modelado_ingresos(p.ingresos) 

# Estimación del punto donde se maximizan los beneficios para el precio óptimo calculado
max_beneficios = modelado_beneficios(p.beneficios) # beneficios with estimated optimum precio

Visualizar el modelo

# Generar la recta donde se optimizan los ingresos
 df_ingresos = data.frame(p.ingresos, max_ingresos)

# Generar la recta donde se optimizan los beneficios
 df_beneficios = data.frame(p.beneficios, max_beneficios)

# Visualizar el modelo
 ggplot(data = data.frame(precio = 0)) +
 stat_function(fun = modelado_ingresos, mapping = aes(x = precio, color = 'Modelado de ingresos')) +
 stat_function(fun = modelado_beneficios, mapping = aes(x = precio, color = 'Modelado de beneficios')) +

geom_segment(aes(x = p.ingresos, y = max_ingresos, xend = p.ingresos, yend = 0), data = df_ingresos) +
 geom_segment(aes(x = p.beneficios, y = max_beneficios, xend = p.beneficios, yend = 0), data = df_beneficios) +

scale_x_continuous(limits = c(4, 12)) +
 labs(title = 'Optimizacion de beneficios e ingresos', colour = '') +
 ylab('resultados')