#! /usr/bin/env python3
#
def ford_model ( ):

#*****************************************************************************80
#
## ford_model() analyzes the Ford data.
#
#  Modified:
#
#    30 January 2022
#
#  Author:
#
#    John Burkardt
#
  import matplotlib.pyplot as plt
  import numpy as np

  print ( "ford_model:" )
#
#  Get the data.
#
  filename = 'ford_data.txt'
  data = np.loadtxt ( filename )
  n, d = np.shape ( data )
  print ( "  Data from " + filename + " involves n =", n, "items with dimension d =", d )
#
#  Analyze the data
#
  print ( "" )
  print ( "  Data statistics:" )
  print ( "    Min      = ", np.min ( data, axis = 0 ) )
  print ( "    Max      = ", np.max ( data, axis = 0 ) )
  print ( "    Range    = ", np.max ( data, axis = 0 ) - np.min ( data, axis = 0 ) )
  print ( "    Mean     = ", np.mean ( data, axis = 0 ) )
  print ( "    Variance = ", np.var ( data, axis = 0 ) )
#
#  Display the mileage and price data.
#
  plt.clf ( )
  plt.plot ( data[:,1], data[:,2], 'ro', markersize = 8 )
  plt.grid ( True )
  plt.xlabel ( '<--- Mileage --->' )
  plt.ylabel ( '<--- Asking price (\$)--->' )
  plt.title ( 'Ford Escort Price as Function of Mileage', fontsize = 16 )
  filename = 'ford_mileage_price.png'
  plt.savefig ( filename )
  print ( '  Graphics saved as "%s"' % ( filename ) )
  plt.show ( )
  plt.clf ( )
#
#  Display the year and price data.
#
  plt.clf ( )
  plt.plot ( data[:,0], data[:,2], 'ro', markersize = 8 )
  plt.grid ( True )
  plt.xlabel ( '<--- Year --->' )
  plt.ylabel ( '<--- Asking price (\$)--->' )
  plt.title ( 'Ford Escort Price as Function of Year', fontsize = 16 )
  filename = 'ford_year_price.png'
  plt.savefig ( filename )
  print ( '  Graphics saved as "%s"' % ( filename ) )
  plt.show ( )
  plt.clf ( )
#
#  Create a column vector of 1's.
#  Create a new data vector X which begins with the column of 1's.
#
  X = np.insert ( data[0:n,0:d-1], 0, np.ones( n ), axis = 1 )
  y = data[:,d-1]
#
#  Normal equations.
#
  XTX = np.matmul ( np.transpose ( X ), X )
  XTy = np.matmul ( np.transpose ( X ), y )
  c1 = np.linalg.solve ( XTX, XTy )
  r1 = np.matmul ( X, c1 ) - y
  mse1 = 1 / n * sum ( r1**2 )
  print ( '  c1=', c1 )
  print ( '  mse1=', mse1 )
#
#  lstsq
#
  c4 = np.linalg.lstsq ( X, y, rcond = None )[0]
  r4 = np.matmul ( X, c4 ) - y
  mse4 = 1 / n * sum ( r4**2 )
  print ( '  c4=', c4 )
  print ( '  mse4=', mse4 )
#
#  Display the mileage and price data.
#
  plt.clf ( )
  plt.plot ( data[:,1], data[:,2], 'ro', markersize = 8 )
  m1 = np.min ( data[:,1] )
  m2 = np.max ( data[:,1] )
  p1 = c1[0] + c1[1] * np.mean ( data[:,0] ) + c1[2] * m1
  p2 = c1[0] + c1[1] * np.mean ( data[:,0] ) + c1[2] * m2
  plt.plot ( [m1,m2], [p1,p2], 'b-', markersize = 8 )
  plt.grid ( True )
  plt.xlabel ( '<--- Mileage --->' )
  plt.ylabel ( '<--- Asking price (\$)--->' )
  plt.title ( 'Ford Escort Price as Function of Mileage', fontsize = 16 )
  filename = 'ford_mileage_price_predicted.png'
  plt.savefig ( filename )
  print ( '  Graphics saved as "%s"' % ( filename ) )
  plt.show ( )
  plt.clf ( )
#
#  Display the year and price data.
#
  plt.clf ( )
  plt.plot ( data[:,0], data[:,2], 'ro', markersize = 8 )
  y1 = np.min ( data[:,0] )
  y2 = np.max ( data[:,0] )
  p1 = c[0] + c[1] * y1 + c[2] * np.mean ( data[:,1] )
  p2 = c[0] + c[1] * y2 + c[2] * np.mean ( data[:,1] )
  plt.plot ( [y1,y2], [p1,p2], 'b-', markersize = 8 )
  plt.grid ( True )
  plt.xlabel ( '<--- Year --->' )
  plt.ylabel ( '<--- Asking price (\$)--->' )
  plt.title ( 'Ford Escort Price as Function of Year', fontsize = 16 )
  filename = 'ford_year_price_predicted.png'
  plt.savefig ( filename )
  print ( '  Graphics saved as "%s"' % ( filename ) )
  plt.show ( )
  plt.clf ( )

  return

if ( __name__ == "__main__" ):
  ford_model ( )