The Mad Economist

F# - Fluent NHibernate - POCO Class Object

2009-10-13T16:44:59Z

The next thing I want to try is to use an F# Class (instead of C#) with Fluent NHibernate.

In my last post I worked with a C# Class & Fluent NHibernate. The C# Class looked as follows:

------------------------------------------------------------------
using System;
using System.Collections.Generic;
using System.Linq;
using System.Text;

namespace StockPrices {

public class COMPANY {

        public virtual int Id { get; set; }
        public virtual string COMPANY_NAME { get; set; }
        public virtual string COMPANY_TICKER { get; set; }
        public virtual int ASSET_TYPE_ID { get; set; }
        public virtual int IS_ACTIVE { get; set; }

}

}
------------------------------------------------------------------

So how do we make the same class in F#?

The first thing I tried was this:

------------------------------------------------------------------
#light

namespace StockPrices

open System
open System.Collections.Generic

// ----------------- 1
type COMPANY = class

        [<DefaultValue>] val mutable Id : int
        [<DefaultValue>] val mutable COMPANY_NAME : string
        [<DefaultValue>] val mutable COMPANY_TICKER : string
        [<DefaultValue>] val mutable ASSET_TYPE_ID : int
        [<DefaultValue>] val mutable IS_ACTIVE : int

        new() = {}

    end
------------------------------------------------------------------

This is the simplest way of setting up a class in F#.

In this example, I use the default constructor -- these are the open parentheses in the New() statement. Therefore, my properties are not initialized. In F# this is a no no and is not allowed. When you do not initialize a variable you need to use the [<DefaultValue>] attribute to initialize the value to Zero, it that is possible.

From Microsoft: "The DefaultValue attribute is required on explicit fields in class types that have a primary constructor. This attribute specifies that the field is initialized to zero. The type of the field must support zero-initialization."

So, the above statement gives me a class with five public variables -- the integers are initialized with Zeros; the strings are initialized with Nulls.

The following is my main F# program. Note that my COMPANY class is in a module called StockPrices that I open at the start of the program.

-----------------------------------------------------------------------------------------------
#light
open System
open System.Collections.Generic
open System.IO
open StockPrices

open FluentNHibernate.Automapping
open FluentNHibernate

let properties = new Dictionary<string, string>()

properties.Add("connection.provider", "NHibernate.Connection.DriverConnectionProvider")
properties.Add("dialect", "NHibernate.Dialect.MsSql2000Dialect")
properties.Add("connection.driver_class", "NHibernate.Driver.SqlClientDriver")
properties.Add("show_sql", "true")

properties.Add("proxyfactory.factory_class", "NHibernate.ByteCode.Castle.ProxyFactoryFactory, NHibernate.ByteCode.Castle")

let connString = "server='BIG_ROCK\LOGGERSEDGE';Initial Catalog=SMDATA;User ID=sa;Password=XXXXX"
properties.Add("connection.connection_string", connString)

//-> The New
let autoMappings = (FluentNHibernate.Automapping.AutoMap.AssemblyOf<StockPrices.COMPANY>()).Where(fun t -> (t.Namespace="StockPrices"))

let aConfig = (new NHibernate.Cfg.Configuration()).AddProperties(properties).AddAutoMappings(autoMappings)

let sessionFactory = aConfig.BuildSessionFactory()

let aSession = sessionFactory.OpenSession()

aSession.BeginTransaction()

//AT&T is Id 100
let coID = 100

let someObj = aSession.Load(typeof<StockPrices.COMPANY>, coID) :?> StockPrices.COMPANY

printfn "Company Name: %s, Ticker: %s" someObj.COMPANY_NAME someObj.COMPANY_TICKER

aSession.Close()

let userresp = Console.ReadLine()
-----------------------------------------------------------------------------------------------

This program craps out on the line:

let aConfig = (new NHibernate.Cfg.Configuration()).AddProperties(properties).AddAutoMappings(autoMappings)

with the error:

{"(XmlDocument)(2,4): XML validation error: The element 'class' in namespace 'urn:nhibernate-mapping-2.2' has incomplete content. List of possible elements expected: 'meta, subselect, cache, synchronize, comment, tuplizer, id, composite-id' in namespace 'urn:nhibernate-mapping-2.2'."}

To be honest, I really have no clue what this all means.

So, I decide to take another tack. I try an F# class that actually uses the 'member' syntax:

-----------------------------------------------------------------------------------------------
    type COMPANY() = class

        let mutable _Id : int = 0
        let mutable _COMPANY_NAME : string = ""
        let mutable _COMPANY_TICKER : string = ""
        let mutable _ASSET_TYPE_ID : int = 0
        let mutable _IS_ACTIVE : int = 0

        member x.Id with get() = _Id and set(v) = _Id <- v
        member x.COMPANY_NAME with get() = _COMPANY_NAME and set(v) = _COMPANY_NAME <- v
        member x.COMPANY_TICKER with get() = _COMPANY_TICKER and set(v) = _COMPANY_TICKER <- v
        member x.ASSET_TYPE_ID with get() = _ASSET_TYPE_ID and set(v) = _ASSET_TYPE_ID <- v
        member x.IS_ACTIVE with get() = _IS_ACTIVE and set(v) = _IS_ACTIVE <- v

end
-----------------------------------------------------------------------------------------------

I rerun the my program.

This time, I make it one more line before it crashes on the line:

let sessionFactory = aConfig.BuildSessionFactory()

This time I get this big ugly:

NHibernate.InvalidProxyTypeException was unhandled
Message="The following types may not be used as proxies:
StockPrices.COMPANY: method get_Id should be 'public/protected virtual' or 'protected internal virtual'
StockPrices.COMPANY: method set_Id should be 'public/protected virtual' or 'protected internal virtual'
StockPrices.COMPANY: method get_COMPANY_NAME should be 'public/protected virtual' or 'protected internal virtual'
StockPrices.COMPANY: method set_COMPANY_NAME should be 'public/protected virtual' or 'protected internal virtual'
StockPrices.COMPANY: method get_COMPANY_TICKER should be 'public/protected virtual' or 'protected internal virtual'
StockPrices.COMPANY: method set_COMPANY_TICKER should be 'public/protected virtual' or 'protected internal virtual'
StockPrices.COMPANY: method get_ASSET_TYPE_ID should be 'public/protected virtual' or 'protected internal virtual'
StockPrices.COMPANY: method set_ASSET_TYPE_ID should be 'public/protected virtual' or 'protected internal virtual'
StockPrices.COMPANY: method get_IS_ACTIVE should be 'public/protected virtual' or 'protected internal virtual'
StockPrices.COMPANY: method set_IS_ACTIVE should be 'public/protected virtual' or 'protected internal virtual'
StockPrices.COMPANY: field _Id should not be public nor internal
StockPrices.COMPANY: field _COMPANY_NAME should not be public nor internal
StockPrices.COMPANY: field _COMPANY_TICKER should not be public nor internal
StockPrices.COMPANY: field _ASSET_TYPE_ID should not be public nor internal
StockPrices.COMPANY: field _IS_ACTIVE should not be public nor internal"

This message provides a few clues. The messages seem to indicate that my class needs to use virtual properties (ala C#). I also assume this to mean that my F# syntax does not yield virtual properties.

So, I need to poke around a find how to create virtual properties in F#

This is what I came up with:

abstract Id : int with get, set
default x.Id with get() = _Id and set(v) = _Id <- v

abstract COMPANY_NAME : string with get, set
default x.COMPANY_NAME with get() = _COMPANY_NAME and set(v) = _COMPANY_NAME <- v

abstract COMPANY_TICKER : string with get, set
default x.COMPANY_TICKER with get() = _COMPANY_TICKER and set(v) = _COMPANY_TICKER <- v

        abstract ASSET_TYPE_ID : int with get, set
        default x.ASSET_TYPE_ID with get() = _ASSET_TYPE_ID and set(v) = _ASSET_TYPE_ID <- v

        abstract IS_ACTIVE : int with get, set
        default x.IS_ACTIVE with get() = _IS_ACTIVE and set(v) = _IS_ACTIVE <- v

end
-----------------------------------------------------------------------------------------------

It is very similar to my class above, however, it uses the keyword 'abstract' to create an abstract property. The 'default' keyword (replacing the member keyword) then heads of the implementation of the property in the current class.

This is from MSDN:

"Properties can be abstract. As with methods, abstract just means that there is a virtual dispatch associated with the property. Abstract properties can be truly abstract, that is, without a definition in the same class. The class that contains such a property is therefore an abstract class. Alternatively, abstract can just mean that a property is virtual, and in that case, a definition must be present in the same class. Note that abstract properties must not be private, and if one accessor is abstract, the other must also be abstract."http://msdn.microsoft.com/en-us/library/dd483467(VS.100).aspx

If one were to refer to a POCO class in F#, this would probably be along the lines you would expect.

OK, let's give this puppy a twirl!

Kaboom!

I crash again on the same line as before.

This is my error:

NHibernate.InvalidProxyTypeException was unhandled
Message="The following types may not be used as proxies:
StockPrices.COMPANY: field _Id should not be public nor internal
StockPrices.COMPANY: field _COMPANY_NAME should not be public nor internal
StockPrices.COMPANY: field _COMPANY_TICKER should not be public nor internal
StockPrices.COMPANY: field _ASSET_TYPE_ID should not be public nor internal
StockPrices.COMPANY: field _IS_ACTIVE should not be public nor internal"

Shorter, but still a problem. At least I have gotten rid of the 'virtual property' errors and I am now left with a bunch of errors that seem to say that my private variable placeholders are illegal. I find this odd, as I do not get the same message with this type of syntax in VB:

-------------------------------------------------------------------------

Option Explicit On

Namespace StockPrices

Public Class COMPANY

        Private _ID As Integer
        Private _CompanyName As String
        Private _CompanyTicker As String

Public Overridable Property Id() As Integer

            Get
                Id = _ID
            End Get

            Set(ByVal value As Integer)
                _ID = value
            End Set

End Property

Public Overridable Property COMPANY_NAME() As String

            Get
                COMPANY_NAME = _CompanyName
            End Get

            Set(ByVal value As String)
                _CompanyName = value
            End Set

End Property

Public Overridable Property COMPANY_TICKER() As String

            Get
                COMPANY_TICKER = _CompanyTicker
            End Get

            Set(ByVal value As String)
                _CompanyTicker = value
            End Set

End Property

End Class

End Namespace

-------------------------------------------------------------------------

The above outline for a class works fine with Fluent NHibernate.

What to do?

Well, I know my class would be OK if I could just get rid of the proxy errors. I mean these variables should be irrelevant, no?

So, I found a way to suppress them:

properties.Add("use_proxy_validator", "false")

I just need to add the above line to my properties collection. This kills the proxy validator.

I know this could have unknown repercussions, so if anyone knows a legitimate way around this problem, I would sure like to know.

Low and behold, this actually works.

I get exactly what I expect!

My output:

-------------------------------------------------------------------------

NHibernate: SELECT company0_.Id as Id0_0_, company0_.IS_ACTIVE as IS2_0_0_, company0_.ASSET_TYPE_ID as ASSET3_0_0_, company0_.COMPANY_TICKER as COMPANY4_0_0_, company0_.COMPANY_NAME as COMPANY5_0_0_ FROM [COMPANY] company0_ WHERE company0_.Id=@p0;@p0 = 100
Company Name: AT&T, Ticker: T

-------------------------------------------------------------------------

We are now good!

F# - Fluent NHibernate RTM

2009-09-22T13:34:09Z

Well, I finally got the courage to download the RTM version of Fluent NHibernate.

I wrote about using Fluent NHibernate in an earlier post (June 28, 2009), so a lot of this is repetitive and some is new.

I have highlighted in yellow the important stuff.

First, I create a simple class

------------------------------------------------------------------
using System;
using System.Collections.Generic;
using System.Linq;
using System.Text;

namespace StockPrices {

public class COMPANY {

}

}
------------------------------------------------------------------

The COMPANY class has only five properties: Id, COMPANY_NAME and COMPANY_TICKER, ASSET_TYPE_ID and IS_ACTIVE

I next create a SQL Server table that is iso-morphic to my class.

CREATE TABLE [COMPANY] (
[Id] [int] IDENTITY (1, 1) NOT NULL ,
[COMPANY_NAME] [varchar] (100) NOT NULL ,
[COMPANY_TICKER] [varchar] (50) NOT NULL ,
[ASSET_TYPE_ID] [int] NOT NULL ,
[IS_ACTIVE] [int] NOT NULL ,
CONSTRAINT [PK_COMPANY] PRIMARY KEY CLUSTERED
(
[Id]
) ON [PRIMARY]
) ON [PRIMARY]
GO

OK, now I can start to work in F#

I now a new project.

I first add a reference to my StockPrices class and also an open statement:

#light
open System
open System.Collections.Generic
open System.IO
open StockPrices

I now need to add references to FluentNHibernate and NHibernate (I do believe that you need both)

I then add two open statements into my code for FluentNHibernate

------------------------------------------------------------------
open FluentNHibernate

//Old
//open FluentNHibernate.AutoMap

//New RTM
open FluentNHibernate.Automapping

------------------------------------------------------------------

Note the change in the namespace "AutoMap" has been changed to "AutoMapping".

The next thing I do is declare a dictionary object to hold a set of configuration properties. I then add a set of attributes.

let properties = new Dictionary<string, string>()

let connString = "server='BIG_ROCK\LOGGERSEDGE';Initial Catalog=SMDATA;User ID=sa;Password=XXXXXX"
properties.Add("connection.connection_string", connString)

Not much to say here, other than I am using SQL Server 2000, so I have to tell nHibernate to use the MsSql2000Dialect. I set my connection string and add it to my set of properties.

In the RTM version, I also had to add the following property:

properties.Add("proxyfactory.factory_class", "NHibernate.ByteCode.Castle.ProxyFactoryFactory, NHibernate.ByteCode.Castle")

I did not need this before, but you need it now. Otherwise you get this big ugly:

The ProxyFactoryFactory was not configured.
Initialize 'proxyfactory.factory_class' property of the session-factory configuration section with one of the available NHibernate.ByteCode providers.
Example:
<property name='proxyfactory.factory_class'>NHibernate.ByteCode.LinFu.ProxyFactoryFactory, NHibernate.ByteCode.LinFu</property>
Example:
<property name='proxyfactory.factory_class'>NHibernate.ByteCode.Castle.ProxyFactoryFactory, NHibernate.ByteCode.Castle</property>

To be honest, I have not explored the guts of NHibernate to figure this all out. I just know I could not get my older version to run without the proxyfactory.factory_class property.

Next, come the autoMappings. This has also changed.

//-> The Old
//let autoMappings = (AutoPersistenceModel.MapEntitiesFromAssemblyOf<StockPrices.COMPANY>()).Where(fun t -> (t.Namespace = "StockPrices"))

//-> The New
let autoMappings = (FluentNHibernate.Automapping.AutoMap.AssemblyOf<StockPrices.COMPANY>()).Where(fun t -> (t.Namespace="StockPrices"))

Not too much difference here -- basically the method "MapEntitiesFromAssemblyOf" has been replaced with "AutoMap.AssemblyOf". It just took a while to figure it all out.

The point of this statement is that I want FluentNHibernate to map my Company class object and my database table automatically, obviating the need for a mapping interface xml file.

This part:

let autoMappings = FluentNHibernate.Automapping.AutoMap.AssemblyOf<StockPrices.COMPANY>()

says map my object StockPrices.COMPANY to the database. I then need to add a where clause that tells FluentNHibernate what namespace to use. You should check out the FluentNHiernate site to get a more technical understanding of what is going on in the background.

I also discovered I need to wrap this:

FluentNHibernate.Automapping.AutoMap.AssemblyOf<StockPrices.COMPANY>()

in parentheses like this:

(FluentNHibernate.Automapping.AutoMap.AssemblyOf<StockPrices.COMPANY>())

in order to get intellisense to work. Also, it generated this convoluted error:

"Successive arguments should be separated by spaces or tupled, and arguments involving function or method applications should be parenthesized"

It took me a while to figure out what this meant. I guess I must be stupid. Anyway, the where clause:

.Where(fun t -> t.Namespace = "StockPrices")

is an anonymous or lambda function. The function returns true if the namespace of an entity equals the name "StockPrices" where my COMPANY class resides.

F# provides a way to define a nameless function using the keyword fun. This type of function receives just one input value and returns just one output value. Generally, if a function is to be passed as an argument to another function (as in the case here), then often you don’t need to give it a name of its own. These functions are referred to as anonymous functions and sometimes called lambda functions or even just lambdas.

The guts of the function above:

t.Namespace = "StockPrices"

could have just as easily have been written:

if t.Namespace = "StockPrices" then true else false

but that would not have been as cool.

My next line:

let aConfig = (new NHibernate.Cfg.Configuration()).AddProperties(properties).AddAutoMappings(autoMappings)

does quite a bit of work: It

(1) Instantiates an NHibernate Configuration;
(2) Add my dictionary lit of properties to the Configuration; and most importantly
(3) Adds my autoMappings

Again, according to Gregory, AddAutoMappings substitutes for AddAssembly (used in regular NHibernate). This allows us to stop NHibernate from looking for hbm.xml files, and use our auto mapped entities instead.

The next block of code opens a NHibernate session and looks up AT&T in my database. It then prints the company name and ticker to the console, giving us the gratification of knowing that something works.

let sessionFactory = aConfig.BuildSessionFactory()

let aSession = sessionFactory.OpenSession()

aSession.BeginTransaction()

//AT&T is Id 100
let coID = 100

let someObj = aSession.Load(typeof<StockPrices.COMPANY>, coID) :?> StockPrices.COMPANY

printfn "Company Name: %s, Ticker: %s" someObj.COMPANY_NAME someObj.COMPANY_TICKER

aSession.Close()

let userresp = Console.ReadLine()

The one weird piece of code is this:

let someObj = aSession.Load(typeof<StockPrices.COMPANY>, coID) :?> StockPrices.COMPANY

should return a type of StockPrices.COMPANY, but it does not. -- it returns a generic Object. I have not figured out why. It certainly does in VB.

So, I have to cast a generic Object into StockPrices.COMPANY so I could actually use it. In F#, there is an operator ":?>" known as on downcast operator the uses the syntax:

x :?> T

According to our buddies at MSFT:

"The :?> operator performs a dynamic cast, which means that the success of the cast is determined at run time. A cast that uses the :?> operator is not checked at compile time; but at run time, an attempt is made to cast to the specified type. If the object is compatible with the target type, the cast succeeds. If the object is not compatible with the target type, the runtime raises an InvalidCastException."

Anyhow, this operator allows me to cast the result into a type StockPrices.COMPANY.

I then print the output to the console. Note: I also get my sql query generated by NHibernate. This is the result of setting the show_sql property to true: properties.Add("show_sql", "true").

Next Task: Use an F# Class instead of a C# Class

Here is the full F# code:

-----------------------------------------------------------------

#light
open System
open System.Collections.Generic
open System.IO
open StockPrices

open FluentNHibernate.Automapping
open FluentNHibernate

let properties = new Dictionary<string, string>()

properties.Add("proxyfactory.factory_class", "NHibernate.ByteCode.Castle.ProxyFactoryFactory, NHibernate.ByteCode.Castle")

let connString = "server='BIG_ROCK\LOGGERSEDGE';Initial Catalog=SMDATA;User ID=sa;Password=XXXXXX"
properties.Add("connection.connection_string", connString)

//-> The Old
//let autoMappings = (AutoPersistenceModel.MapEntitiesFromAssemblyOf<StockPrices.COMPANY>()).Where(fun t -> (t.Namespace = "StockPrices"))

//-> The New
let autoMappings = (FluentNHibernate.Automapping.AutoMap.AssemblyOf<StockPrices.COMPANY>()).Where(fun t -> (t.Namespace="StockPrices"))

let aConfig = (new NHibernate.Cfg.Configuration()).AddProperties(properties).AddAutoMappings(autoMappings)

let sessionFactory = aConfig.BuildSessionFactory()

let aSession = sessionFactory.OpenSession()

aSession.BeginTransaction()

//AT&T is Id 100
let coID = 100

let someObj = aSession.Load(typeof<StockPrices.COMPANY>, coID) :?> StockPrices.COMPANY

printfn "Company Name: %s, Ticker: %s" someObj.COMPANY_NAME someObj.COMPANY_TICKER

aSession.Close()

let userresp = Console.ReadLine()

-----------------------------------------------------------------

Excuse Me?

2009-09-16T20:00:10Z

Well, now I have definite evidence that there is no intelligent life on earth.

Last week an appeals court threw out Lucent's $358 million award against Microsoft for infringing on a patent allowing MS Outlook users to enter a date by clicking on a calendar. Go here http://www.channelinsider.com/c/a/Errata/Microsoft-358-Million-Damage-Award-Overturned-620309?.

Does anyone else think this is galactically stupid?

I have have just three questions.

(1) What idiots on a jury think clicking on a calendar is worth $358 Million?

(2) Does Lucent think that because they're such a failure they have to resort to frivolous lawsuits to make money?

(3) What on earth is the patent office thinking by granting a patent for clicking? Next will come addition or multiplication.

It seems like the only one with any intelligence was the judge!

Come on you guys. Fire all your lawyers, hire engineers and get working on the cool stuff.

F# - Classes (In progress)

2009-09-05T10:07:02Z

So, today I decided I would try to figure out how classes work in F#.

I have decided to make a Bank object. The first thing I tried is this:

#light

type Bank() = class
nCustomers : int
nEmployees : int
end

This does not work! I get the error:

Program.fs(6,5): error FS0010: Unexpected identifier in member definition.

Who knows what that means.

So, I start poking around and find that I need to use the 'val' key word to create a public property. So the next thing I tried is this:

type Bank() = class
val nCustomers : int
val nEmployees : int
end

This also does not work! I get the error:

Program.fs(6,9): error FS0191: Uninitialized 'val' fields in implicit construction types must be mutable and marked with the '[<DefaultValue>]' attribute. Consider using a 'let' binding instead of a 'val' field

After poking around on Mister Softy, I find I need to do it like this:

type Bank() = class
[<DefaultValue>] val mutable nCustomers : int
[<DefaultValue>] val mutable private nEmployees : int
end

So this is why:

(1) I need the val keyword in conjunction with the mutable keyword in order to make the property public. That is the default. Note that my property nEmployees uses the keyword 'private' to make the property private. The let construct will not work if you want the property to be public. Let bindings are always private.

(2) I really do need both the val and mutable keywords.

(3) In this example, I use the default constructor -- these are the open parentheses in the Bank() statement. Therefore, my properties are initialized. In F# this is a no no and is not allowed. When you do not initialize a variable you need to use the Default attribute to initialize the value to Zero, it that is possible. From Microsoft:

"The DefaultValue attribute is required on explicit fields in class types that have a primary constructor. This attribute specifies that the field is initialized to zero. The type of the field must support zero-initialization."

So, the above statement gives me a class with one public variable (nCustomers) initialized at zero and one private variable, also initialized at zero.

The following code snippet shows how I instantiate a Bank and assign it a nCustomers:

let BofA = new Bank()

//The following does not work
//let BofA.nCustomers = 4300

BofA.nCustomers <- 4300

Now what if I do want to use a constructor instead of the default constructor. I can create a constructor thusly:

type Bank = class

[<DefaultValue>] val mutable nCustomers : int
[<DefaultValue>] val mutable private nEmployees : int

new(nCust) = {nCustomers = nCust}

end

In the above snippet, I now create a constructor and use it to assign one of my properties nCustomers. The above, however, does not work. I get this error:

Program.fs(10,18): error FS0191: Extraneous fields have been given values

I guess F# can't deal with both a constructor and a default initialization at the same time. So, I need to write it as such:

type Bank = class

val mutable nCustomers : int
[<DefaultValue>] val mutable private nEmployees : int

new(nCust) = {nCustomers = nCust}

end

I need to take out the attribute [<DefaultValue>] if the property is assigned by the constructor. I guess that makes sense.

So, now I can instantiate my Bank and then cut the numbers of customers as as such:

let BofA = new Bank(5400)

BofA.nCustomers <- 4300

If I want to initialize both my properties, I would make by constructor as follows:

type Bank = class

val mutable nCustomers : int
val mutable private nEmployees : int

new(nCust, nEmp) = {nCustomers = nCust; nEmployees = nEmp}

end

You separate your initialization statements with a semi-colon.

You can also set up a class using an implicit constructor as such:

type Bank(nCust : int, nEmp : int) = class

    let mutable _Cust : int = nCust
    let mutable _Emp : int = nEmp

    member x.nCustomers with get() = _Cust and set(v) = _Cust <- v
    member x.nEmployees with get() = _Emp and set(v) = _Emp <- v

end

This format allows you to put the constructor in the class header. You then set up get and set statements as needed for your class properties. Note that you don't need both a get and set statement. If the property is read only, you would use the syntax:

member x.nCustomers with get() = _Cust

Or, if you limited the property to write only the syntax would be:

member x.nCustomers with set(v) = _Cust <- v

So what is the "x." doodad for? Well, I can tell you one thing, if you don't use it you will get an invalid syntax error.

The x. really is the class "self" identifier, like Me or this. Microsoft points out that is great because "Unlike in other .NET languages, you can name the self identifier however you want; you are not restricted to names such as self, Me, or this."

Inheritance

Inheritance works in F#. Using an explicit constructor we can create a base class (e.g. Bank) and then have another class, say CitiBank inherit the properties and methods of the Bank class:

type Bank = class

val mutable nCustomers : int
val mutable nEmployees : int

new(nCust, nEmp) = {nCustomers = nCust; nEmployees = nEmp}

end

type CitiBank = class
inherit Bank

val nDeposits : float

    new(nCust : int, nEmp : int, xDep : float) =
        {inherit Bank(nCust, nEmp);
         nDeposits = xDep}

end

You just need to use the key word inherit (lower case) in the class declaration. We can then add additional properties (e.g. deposits) to our new class. Note the syntax of the new constructor for the CitiBank class -- you basically call the constructor for the base class, passing along the arguments from the constructor of the higher level class.

If you want to use an implicit constructor, this syntax should work:

type Bank(nCust : int, nEmp : int) = class

    let mutable _Cust : int = nCust
    let mutable _Emp : int = nEmp

    member x.nCustomers with get() = _Cust and set(v) = _Cust <- v
    member x.nEmployees with get() = _Emp and set(v) = _Emp <- v

end

type FedBank(nCust : int, nEmp : int, xDep : float) = class
inherit Bank(nCust, nEmp)

let mutable _Dep : float = xDep
member x.nDeposits with get() = _Dep and set(v) = _Dep <- v

end

Again, you call the constructor for the base class in the constructor for the higher level class.

Adding a Method

So, now let's add a method called 'AddDeposit'. We use the member keyword and then define our function. The function takes a float and adds it to the existing deposit balance. The function doesn't return anything.

type FedBank(nCust : int, nEmp : int, xDep : float) = class
inherit Bank(nCust, nEmp)

let mutable _Dep : float = xDep
member x.nDeposits with get() = _Dep and set(v) = _Dep <- v

    member x.AddDeposit(aDeposit : float) =
        _Dep <- _Dep + aDeposit
        ()

end

I need to use the 'mutable' syntax because the deposit balance will be constantly changing.

We can now see if it works:

let BofA = new FedBank(5400, 56, 32212.21)

let dep0 = BofA.nDeposits

BofA.AddDeposit(5000.0)

let dep1 = BofA.nDeposits

This does work.

Note that if I had wanted to return the current deposit balance, I would have written my function:

    member x.AddDeposit(aDeposit : float) =
        _Dep <- _Dep + aDeposit
        _Dep

Railroad Traffic Data - WE 8/29/2009

2009-09-03T11:07:14Z

Well, I just got a look at the rail traffic data for the week ending 8/29/2009.

We actually see a fairly sizable increase -- all non-intermodal traffic increased by 2.1% over the previous week:

Two percent may not sound like a lot, but in the world of railroads it is a fairly big deal. Total carloads are now 285,580, well off their low of around 250,000 just a few weeks ago.

Let's look at coal:

Coal loads were up just about the same 2.1 percent. This stands to reason as coal makes up nearly half of all carloads.

Chemicals also showed a fairly good sized jump after being fairly flat:

The week-over-week increase for chemicals was 6.25 percent. Carloads reached their highest level in nearly a year.

Intermodal traffic also showed a fairly good jump, after being virtually flat. Intermodal traffic increased by 4.84 percent -- again a sizable increase for just one week.

These are all pretty good numbers. So, what can we say? This week's numbers provide further assurance that things are starting to pick up.

Where is the Economy - Railroad Traffic 8/22/09

2009-08-28T07:22:26Z

Here is the chart for railroad traffic for the week ended 8/22/09

The data are for all non-intermodal rail carloads. I have been following this statistic as an indicator of the state of the economy. It appears that the trend started about 5 or 6 weeks ago of increasing rail traffic is holding to form. Total non-intermodal traffic stands at about 280,000 carloads versus 250,000 carloads (an increase of about 12 percent). While I was concerned a few weeks ago that the uptick was a blip, I am more confident now that we are seeing is the real deal. I should note that traffic this time last year was on order of 325,000 carloads, so we are far from being out of the woods.

Intermodal carloads are shown below:

The pattern here is not as positive. Though there has been an uptick in the past few weeks, the most recent observations are flat.

The biggest user of the railroads is the coal industry. The recent trend in coal shipments is shown below:

The trend is positive, but still slow. Like a big oil tanker, the US economy is large, unwieldy, and hard to turn. But, methinks, once the momentum has turned the progress will be tough to reverse.

Cash for clunkers?

The following chart show carloads for motor vehicles and equipment.

The recent trend does seem positive, but I don't know if that is a normal seasonal pattern or some offshoot of the cash for clunkers program.

Sugar Shortage #2

2009-08-23T10:28:48Z

In my last post, I complained about the high cost of sugar quotas to US consumers, measuring about $1.8 billion a year. But, what does this have to do with a sugar shortage? In the last few months, the world sugar market has changed dramatically. Although historically the world price of sugar has been far below the US price, the gap has closed considerably. Below, I show the same price chart from my prior post:

Currently, the world price is 21.84 cents/pound and the US price is 26.45 cents/pound. The world price is actually above the target market price set by farm legislation and implemented by the USDA. Ironically, the cost of the subsidy paid to sugar producers is dramatically lower in the current year. The high world price indicates that the supply/demand relationship in the world market has shifted such that its is now a sellers market -- there is not enough sugar on the world market to supply consuming countries.

According to data from the USDA, world production of sugar declined from 166.5 million tons (metric) in the 2007/8 fiscal year to 148.7 million tons (metric) in fiscal 2008/9, a decline of roughly 10 percent. This decline has been matched by an increase in world consumption by about 1.5 percent. What was a market that had bountiful supplies in 2007/8 now has a very tight supply/demand situation.

So, where is sugar produced:

The above chart shows sugar production for the 2008/9 fiscal year. Brazil is by far the largest produce of sugar (the above is for sugar alone, exclusive of cane used for ethanol). India is the second largest producer, followed closely by China. It also turns out that India is also the biggest consumer of sugar (in total) ahead of both China and the United States.

Now the above charts shows absolute sugar production in 2008/9. The following shows the year-to-year change between 2007/8 and 2008/9.

It turns out that the bulk in the world decline in sugar production comes from a really bad crop in India. Of the global reduction in sugar production of 17.7 million tons (metric), India accounts for 11.85 million tons (tons). In fact India has switched from being a net exporter (5.8 million tons in 2007/8) to being a net importer of sugar. The ramifications of this reduction are apparent in the world sugar market -- prices have risen dramatically.

How does this affect the United States?

In 2007/08, Total US Sugar Supply (exclusive of existing stocks) was 10.772 million tons. Of this amount 76% was supplied domestically, 13% from Tariff Rate Quota (TRQ) imports, 6% from Mexico and 5% from other sources. Tariff Rate Quota (TRQ) imports are those limited by Farm Legislation and implemented by the USDA. These imports amounted to 1.354 million tons -- basically the amount of the limit.

In 2008/09, US sugar producers had a bad year and their production declined by 7% -- about 0.581 million tons. So, right from the get go, the supply/demand situation would be tight anyway. Total supply was about the same in 2008/09 -- 10.762 million tons. So what made up the shortfall?

In 2008/09 TRQ imports were basically flat. Mexico basically drastically increased there imports into the US -- growing from 0.694 million tons to 1.450 million tons -- more than double. Mexican sugar is not subject to an import quota and can imported under NAFTA duty free. Methinks that if the red state farmers had had their way in limiting Mexican sugar imports, the pricing situation for sugar would be much worse. The result of the Mexican imports is that the price of sugar in the US has remained fairly flat throughout the year, even though world sugar prices have escalated substantially.

Crunch Time

Right now, the markets are intensely focused on the 2009/10 supply in the US. Domestic production is expected to rise somewhat, basically back to 2007/08 levels. But, overall supply is expected to fall by 0.730 million tons, or about 6.8 percent. How can that be?

The first problem is that TRQ imports are expected to fall by about 0.250 million tons. How can that be. The US limits sugar imports doesn't it? Well the US does limit imports. But it does so on a country-by-country basis. So, if a country like Zimbabwe does not meet its allowable quota of 12,636 metric tons, the quota 'shortfall' just goes into space and is unused. Thus, for example, Brazil has net its quota by the end of June (quotas are set on a fiscal year -- ending September 30th -- basis), it cannot import anymore into the US. The country-by-country limits always mean that there will be an aggregate shortfall. That is, counties that have extra sugar to export can't make up up for those with a short.fall. The result -- depending on the year TRQ imports are always less than what is allowed in aggregate.

So, while the April, 2009 projection of sugar imports for 2009/10 was set at 1.496 million tons, it has been systematically been reduced to 1.182 million tons (a decrease of about 20%) reflecting not a tighter quota limit, but a greater shortfall. Why? Two reasons come to mind. (1) Bad weather & bad harvests or (2) a willingness on producers to sell their goods elsewhere. I suspect it it the latter. Because the world price is basically the same as the US price (as faced by an exporter), why would anyone prefer to sell to the US when they can get the same price from someone else? I mean put yourself in the shoes of some guy in (pick a country) who has been tormented by quotas & paperwork for years. Why would he (or she) voluntarily choose to trade with the US?

Whatever the the rationale, TRQ imports are expected to decrease by 0.250 million tons.

What about Mexico? Mexico is interesting because the country is now producing less sugar than it consumes. So, how did it significantly raise its imports into the US in 2008/09. Answer: It all came out of existing stocks. Mexico has basically sold its entire preexisting (as of the start of 2008/09) stock of sugar to the US. I'm sure they did this because they got a good price. But, the fact of the matter is the cookie jar is empty. Imports from Mexico into the US are expected to fall by 1.285 million tons. Folks, that is a lot of sugar. Mexico will consume what it produces, but there is nothing extra to export.

The USDA also expects that part of the 2009/10 sugar shortfall will be filled by a draw down in US stocks from 1.252 million tons to 0.709 million tons for an extra 0.543 million tons.

These factors all combine to make for a very tight market in the forthcoming fiscal year. This, then, is the focus of the food companies in the letter to the Secretary of Agriculture requesting an increase in the TRQ country limits (see my prior post).

The current US supply/demand balance now reflects the world's supply/demand balance. Currently, the TRQ import quotas are basically non binding and the US must compete for sugar on the world market. If this is in fact the case, it is not immediately obvious that by raising the TRQ quotas (as requested by the food companies), there will be any additional sugar provided into the US market.

A digression on the TRQ program

The TRQ program limits the sugar imports by allocating a quota to 40 countries. The actual quota assigned to a country is based on the total allowable quota and a percentage share assigned to the country. The percentage shares are based on trade data from 1975-81 when the sugar trade was relatively unrestricted.

The chart below shows the 2008/09 TRQ limits for each country compared to their actual imports through July, 2009.

The chart above shows that most of the big producers/importers are nearly maxed out on the quota, with the exception of Argentina which still has a lot of quota room left. There are also a number of smaller producers that have not imported any sugar into the US this year. Overall, after 10 months of the fiscal year, total TRQ imports are running at 71% of quota. By comparison, by the same point in time last fiscal year, imports were running 74% of quota.

Given the above distribution of imports, it appears that for many countries raising the import quotas would have no affect on the sugar supplied to the US market -- many countries won't hit the current quotas. On the other hand, for the big producers like Brazil, The Philippines, Australia, and The Dominican Republic, raising their quotas may indeed yield greater imports because they are likely to max out under their current quotas.

Sugar Quotas #1

2009-08-21T18:10:08Z

I came across a story in the WSJ last week (WSJ August 13th) in which some major food companies were warning of a sugar shortage in the US. Excuse me. A sugar shortage? I'm not making this up. Here's the text of the letter:

August 5, 2009
The Honorable Thomas J. Vilsack
Secretary
U.S. Department of Agriculture
Jamie L. Whitten Federal Building
1400 Independence Avenue, S.W.
Washington, DC 20250

Dear Mr. Secretary:

The organizations and companies below urge you to increase the sugar import quota
immediately. Your experts forecast unprecedented shortages without prompt action.
According to USDA's World Agricultural Supply and Demand Estimates, the United
States will end the next fiscal year with less than 13 days' worth of sugar on hand, unless
imports are increased. If this forecast is accurate, our nation will virtually run out of sugar.

The shortage does not have to happen. The only reason markets are forecast to be
so tight is the restrictive U.S. policy on sugar imports. Imports are subject to restrictive
quotas.

But you have the authority to increase the sugar import quota, and we urge you to do
so immediately, both for the current fiscal year - where high prices already indicate a
painfully tight market - and for the upcoming year.

Without a quota increase, consumers will pay higher prices, food manufacturing jobs
will be at risk and trading patterns will be distorted. Please act now in the interest of all
Americans.

Sincerely,

American Bakers Association
American Beverage Association
Blommer Chocolate
Competitive Enterprise Institute
ConAgra Foods, Inc.
Consumer Federation of America
Council for Citizens Against Government Waste
Emergency Committee for American Trade
General Mills, Inc.
Gonella Frozen Products
Grocery Manufacturers Association
The Hershey Company
Independent Bakers Association
International Dairy Foods Association
Kraft Foods
Krispy Kreme Doughnut Corp.
Mars, Incorporated
National Confectioners Association
Nestle USA
Pepperidge Farm
Sweetener Users Association
Unilever United States, Inc.
U.S. Chamber of Commerce

I say to myself this is bizarre. So, I've spent the last fews days trying to figure out what is going on.

In my sleuthing around, I find out that the U.S. uses price supports, domestic marketing allotments, and tariff-rate quotas (TRQs) to influence the amount of sugar available to the U.S. market. The program is used to support U.S. sugar prices above comparable levels in the world market. In order to keep the program off the Federal Budget, the program relies on import quotas to limit imports, so that a certain target price can be maintained.

The rates for raw cane sugar are:

18 cents per pound in FY 2009,
18.25 cents per pound in FY 2010,
18.50 cents per pound in FY 2011, and
18.75 cents per pound in FY 2012-13.

The rates for refined beet sugar are:

22.9 cents per pound in FY 2009 and
128.5 percent of the loan rate for raw cane sugar in FY 2010-13.

So, I ask myself, what is the world price? In turns out that over the last 15 or so years, the average difference between the world price and the domestic price is about ten cents a pound:

I mean this is about double the world price.

If we look at US Sugar Consumption over the same time period:

It turns out that the US consumes about 900,000 tons a month -- that is a lot of donuts and cokes. You figure in a year that is about 10 million tons. If you figure the US population is 300 million, that means people eat about 67 pounds of sugar a year. Yikes.

So, how much is the subsidy we pay to sugar farmers. I looked at the price difference per month and multiplied it by the monthly US production of sugar:

So, over the last 15 or so years, we as consumers pay sugar producers a subsidy of about $150 Million per month or about $1.8 Billion a year.

So who gets all this dough.

Let's look at who grows sugar cane and sugar beets:

Farm Count, by State:

State	Sugar Beets	Sugar Cane
Caifornia	155	0
Colorado	226	0
Delaware	3	0
Florida	0	108
Georgia	35	0
Hawaii	0	9
Idaho	507	0
Louisiana	0	461
Michigan	737	0
Minnesota	1247	0
Montana	220	0
Nebraska	162	0
North Dakota	553	0
Oregon	73	0
Texas	0	114
Washington	3	0
Wyoming	139	0
Total	4022	692

Sugar Production ('000 Tons), by State:

State	Sugar Beets	Sugar Cane
Caifornia	231	0
Colorado	115	0
Delaware	3	0
Florida	0	1,800
Georgia	*	0
Idaho	845	0
Hawaii	0	200
Louisiana	0	1,400
Michigan	530	0
Minnesota	1,715	0
Montana	176	0
Nebraska	160	0
North Dakota	854	0
Oregon	52	0
Texas	0	152
Washington	12	0
Wyoming	101	0
Total	4,791	3,552

Now let's look at the implied sugar subsidy by multiplying the tons by $0.10/Lb.

Sugar Subsidy ($million), by State:

State	Sugar Beets	Sugar Cane
Caifornia	46.2	0
Colorado	22.9	0
Delaware	*	0
Florida	0	360.0
Georgia	*	0
Idaho	168.9	0
Hawaii	0	40.0
Louisiana	0	280.0
Michigan	106.0	0
Minnesota	343.0	0
Montana	35.2	0
Nebraska	31.9	0
North Dakota	170.9	0
Oregon	10.4	0
Texas	0	30.4
Washington	2.4	0
Wyoming	20.3	0
Total	$958.1	$710.4

If we put this on a 'Per Farm' basis, we can estimate how much subsidy goes to each farm, on average, by state:

Sugar Subsidy ($thound) per Farm, by State:

State	Sugar Beets	Sugar Cane
Caifornia	297.9	0
Colorado	101.5	0
Delaware	*	0
Florida	0	3,333.3
Georgia	*	0
Idaho	333.2	0
Hawaii	0	4,444.4
Louisiana	0	607.4
Michigan	143.8	0
Minnesota	275.0	0
Montana	160.1	0
Nebraska	197.1	0
North Dakota	309.0	0
Oregon	142.3	0
Texas	0	266.7
Washington	802.1	0
Wyoming	145.8	0
Average	238.2	1026.6

Good grief. When all is said and done we are paying sugar beet farmers an average of $238,200/year to grow their crop. We are paying sugar cane farmers average of $1,026,600/year to grow their crop. That's one million bucks a year folks. That's real money. We're talking a serious case of Welfare for Farmers. I need to get in on this racket.

How did I get off on this tangent anyway? I was taliking about a sugar shortage not farm welfare. This post is already too long, that will just have to in the next one.

F# Serial Correlation (2)

2009-08-14T17:34:56Z

In my last post I set up an example dataset and then showed how use the Durbin-Watson statistic to test for serial correlation. Once we have recognized that serial correlation exists, we need to make the appropriate adjustments because OLS is no longer the best estimation technique.

We need to rely on GLS to perform our estimation. The key to this is to identify what the matrix Omega (Z) looks like. In the case where we have first order autocorrelation AR(1), the Omega (Z) matrix looks like:

where the p is first order auto-correlation coefficient -- meaning the correlation between sequential error terms is p (rho). I would use Greek letters but I cant.)

So, we can easily apply GLS, we just need to figure out what p (rho) is.

One way of estimating p (rho) is to use the Durban-Watson statistic. It turns out that the DW Statistic is approximately equal to:

2 * (1 - p). In our regression, we computed a DW statistic of 1.09, which would imply a p (rho) of 0.455.

But, there is another, more interesting method, called the Hildreth-Lu Procedure.

This approach uses a grid search to estimate the most 'likely' value for p (rho). In practice this might involve running 100 regressions, each using a different value of p (rho), such as 0.01, 0.02, 0.03, ...0.99. We would then select the regression that yielded the lowest sum of squared errors. Once we have narrowed down the range, we could run another 100 regressions to further refine our estimate. What a perfect task for F# where we can send off multiple regressions at once and make full use of our dual core or multi-core computers.

We first transform our original data to look like this:

Y* = X* + v

Where

Y* = Y_t - pY_t-1

Each entry in X*: X_kt = X_kt - pX_kt-1

v = e_t - pe_t-1

We then run a regular OLS regression and calculate the sum of squared errors. What we need is a function that will take a value of p (rho), do the data transformation and then run the OLS regression. The function would then return of the sum of squared errors.

Here's the code:

type rssType = {rho : float; rss : float}

let HL_iteration (X0 : float [,] ) (y : float [,] ) (p : float) =

let n1 = n-1

let X2 = Array2D.zeroCreate<float> n1 k
let y2 = Array2D.zeroCreate<float> n1 1

    for i = 0 to (n1-1) do
        y2.[i, 0] <- y.[i+1, 0] - p * y.[i, 0]
        X2.[i, 0] <- 1.0 - p
        X2.[i, 1] <- X0.[i+1, 1] - p * X0.[i, 1]
        X2.[i, 2] <- X0.[i+1, 2] - p * X0.[i, 2]
        X2.[i, 3] <- X0.[i+1, 3] - p * X0.[i, 3]
        X2.[i, 4] <- X0.[i+1, 4] - p * X0.[i, 4]

    let ret, RSS = regress X2 y2
    let retVal = {rho=p; rss=RSS}
    retVal

I first create a type rssType that can hold both rho and the sum of squared errors. The function takes our X and y matrices as input, as well a value for p (rho). We then get the count of observations and reduce it by one because we lose the first observation in the process (the first observation has no lagged value. We then create two arrays, X2 and y2 that will hold our transformed data. We loop over the original data and do the transformation element-by-element. Once we have our transformed data, we pass it to our regular OLS regression method.

We return the values of p (rho) and the sum residual sum of squares (sum or squared errors). The idea here is that we can call this function multiple times for different p's (rho's) and get back the sum of squares. We can then pick the rho with the smallest sum of squares.

This is how we do that:

let HL_iterate (X : float [,] ) (y : float [,] ) (i0 : float) (i1 : float) =

    let incr = (i1 - i0) / 10.0
    let rssList = new List<rssType>()

let asyncList = new List<Async<rssType>>()

    for i = 0 to 9 do
        let (iret : Async<rssType>) = async {return HL_iteration X y ((i0 + ((float) i * incr)))}
        asyncList.Add iret

    let Result =
        async { let! asyncWorkflow = Async.Parallel asyncList
                return asyncWorkflow }

    let rssList = Async.RunSynchronously Result

    for aRSS : rssType in rssList do
        printfn "RSS, rho: %10.5f %10.6E " aRSS.rss aRSS.rho

rssList

This function takes our original X and y matrices, and two end points i0 and i1. The plan here is to take the two end points for rho, say 0.01 and 0.99, and split them into 10 equally spaced values for rho. We then run a regression for each of the ten values and obtain the results. By narrowing then end points we can fine tune our grid search to obtain a more precise estimate of rho. The variable incr is just the range i1 - i0 divided by 10. I then create a list rssList of type rssType to hold my results.

Now, I want to send off all my ten regressions at once, and take advantage of F#'s handling of multi-threaded operations. So, I set up a list of function calls I want to run asynchronously. These will be held in my variable asyncList. I next create a loop with 10 iterations than creates the function call for each of my spaced out values of rho. The value of rho: (i0 + ((float) i * incr) is just the starting point plus i times the value of the i of the increment variable. I add each function call to my asyncList.

I then set up the workflow to process the function calls in parallel and only return when they are all complete (this is the magic of the let! statement).

    let Result =
        async { let! asyncWorkflow = Async.Parallel asyncList
                return asyncWorkflow }

I return the values in the variable asyncWorkflow which is a list of rssType values. Note: this statement does not actually run the functions (regressions). You actually need the following statement:

let rssList = Async.RunSynchronously Result

to run the regressions and fill up the rssList of results.

I then create a main function to call my HL_iterate function:

let rssList = HL_iterate X y 0.01 0.99

for aRSS in rssList do
printfn "RSS, rho: %10.6E %10.5f " aRSS.rss aRSS.rho

let rssArr = Array.sortWith (fun (x : rssType) (y : rssType) -> if y.rss = x.rss then 0 elif y.rss < x.rss then 1 else -1 ) rssList

printfn "Lowest RSS, rho: %10.6E %10.5f " (rssArr.[0].rss) (rssArr.[0].rho)

---------------------------------------------------------------------------

The variable rssList has the results from the ten regressions. I then print them out. My next step is to find the regression with the lowest sum of squares. I do this by invoking the Array method sortWith. This method allows me to sort an array using an anonymous function that I provide.

(fun (x : rssType) (y : rssType) -> if y.rss = x.rss then 0 elif y.rss < x.rss then 1 else -1 )

My function takes the value of rss (one part of the rssType structure) and compares it to the next value. My comparison uses the '<' operator so that I can get a resulting array that is sorted in ascending order. Thus, it is easy to obtain the regression with the lowest sum of squares -- it ends up at the head of the list.

So, let's run it and see what we get:

RSS, rho: 5.018932E+009    0.01000
RSS, rho: 4.636517E+009    0.10800
RSS, rho: 4.348996E+009    0.20600
RSS, rho: 4.156232E+009    0.30400
RSS, rho: 4.058028E+009    0.40200
RSS, rho: 4.054117E+009    0.50000
RSS, rho: 4.144200E+009    0.59800
RSS, rho: 4.328156E+009    0.69600
RSS, rho: 4.606552E+009    0.79400
RSS, rho: 4.981314E+009    0.89200
Lowest RSS, rho: 4.054117E+009    0.50000

I just love to see the both processor cores spike! We see that the lowest rho is in the middle at 0.50. We can fine tune this by setting our end points to 0.5 +/- 0.1 and calling the HL_iterate function again:

HL_iterate X y (rssList.[0].rho - 0.1) (rssList.[0].rho+0.1)

let rssList2 = HL_iterate X y (rssArr.[0].rho - 0.1) (rssArr.[0].rho + 0.1)

for aRSS in rssList2 do
printfn "RSS, rho: %10.6E %10.5f " aRSS.rss aRSS.rho

let rssArr2 = Array.sortWith (fun (x : rssType) (y : rssType) -> if y.rss = x.rss then 0 elif y.rss < x.rss then 1 else -1 ) rssList2

printfn "Lowest RSS, rho: %10.6E %10.5f " (rssArr2.[0].rss) (rssArr2.[0].rho)

We now get our second round of results

RSS, rho: 4.059089E+009    0.40000
RSS, rho: 4.050250E+009    0.42000
RSS, rho: 4.045336E+009    0.44000
RSS, rho: 4.044344E+009    0.46000
RSS, rho: 4.047272E+009    0.48000
RSS, rho: 4.054117E+009    0.50000
RSS, rho: 4.064876E+009    0.52000
RSS, rho: 4.079547E+009    0.54000
RSS, rho: 4.098129E+009    0.56000
RSS, rho: 4.120619E+009    0.58000
Lowest RSS, rho: 4.044344E+009    0.46000

We now add a third iteration, where we adjust the endpoints to be 0.46 +/- 0.01:

let rssList3 = HL_iterate X y (rssArr2.[0].rho - 0.01) (rssArr2.[0].rho + 0.01)

for aRSS in rssList3 do
printfn "RSS, rho: %10.6E %10.5f " aRSS.rss aRSS.rho

let rssArr3 = Array.sortWith (fun (x : rssType) (y : rssType) -> if y.rss = x.rss then 0 elif y.rss < x.rss then 1 else -1 ) rssList3

printfn "Lowest RSS, rho: %10.6E %10.5f " (rssArr3.[0].rss) (rssArr3.[0].rho)

We now get our third round of results:

RSS, rho: 4.044350E+009    0.45000
RSS, rho: 4.044270E+009    0.45200
RSS, rho: 4.044230E+009    0.45400
RSS, rho: 4.044229E+009    0.45600
RSS, rho: 4.044267E+009    0.45800
RSS, rho: 4.044344E+009    0.46000
RSS, rho: 4.044460E+009    0.46200
RSS, rho: 4.044616E+009    0.46400
RSS, rho: 4.044811E+009    0.46600
RSS, rho: 4.045045E+009    0.46800
Lowest RSS, rho: 4.044229E+009    0.45600

So, I think I'll stop here. My best guess as rho is 0.456.

Now I'll just build my Z matrix (shown at the top) and apply my GLS estimation method.

This two-step procedure is commonly referred to as the Prais-Winsten procedure.

It turns out, given the matrix of the form:

that the inverse is easily computed. The inverse of Z is:

Z^-1 = (1 / (1 - p²) * [ 1    -p     0     0     0     0]
                      |-p  1+p²    -p     0     0     0|
                      | 0    -p  1+p²    -p     0     0|
                      | ..   ..    ..    ..    ..    ..|
                      | 0     0     0    -P  1+p²    -p|
                      [ 0     0     0     0   -p     1]

So for our program all we have to do is create Z^-1 based on our estimate of p (rho) and then run our GLS regression

We now add the following to our program:

let rho = rssArr3.[0].rho

//----------------------------------------------------------------------
// Now Build ZI Matrix
//----------------------------------------------------------------------

let dd = (1.0 / (1.0 - Math.Pow(rho,2.0)))

let ZI = Array2D.init<float> n n (fun i j -> if i=j then dd * (1.0 + Math.Pow(rho,2.0)) elif (i=j+1 or j=i+1) then -1.0*rho*dd else 0.0)

ZI.[0,0] <- dd * 1.0
ZI.[n-1,n-1] <- dd * 1.0

let ret = regressGLS X y ZI

We first grab rho estimate from our last Hildreth-Liu iteration. Recall it should be 0.456. We then build a ZI matrix, which is Omega Inverse.

The variable dd is just the multiplier out in front of the matrix (= 1 / (1 - p²)).

ZI is an n x n matrix where the diagonal elements are 1 + p² and the adjacent entries are merely -p. The diagonal corners are set to 1.0. All these elements are multiplied by dd.

I then call my regressGLS method, passing to it our data matrices and Omega inverse (ZI). The full code of RegressGLS is shown below:

--------------------------------------------------------------------------------

let regressGLS (x0 : float [,] ) (y : float [,] ) (ZI : float [,] ) =

let xp = transpose2 x0

    //--------------------------------
    // Now create X'ZI X
    //--------------------------------

let xpx = matMult2 (matMult2 xp ZI) x0

let X = xpx
let XI, ret = matInverse X

    let k = X.GetUpperBound 0
    let ni = 1 + x0.GetUpperBound 0

    for i = 0 to k do
         //printfn "XI Matrix: %10.5f %10.5f " XI.[i,0] XI.[i,1]
         printfn "XI Matrix: %15.10f %15.10f %15.10f %15.10f %15.10f " XI.[i,0] XI.[i,1] XI.[i,2] XI.[i,3] XI.[i,4]

    //-----------------------------------------
    //Now Test
    //Multiply XI and X to calculate Identity Matrix
    //-----------------------------------------
    let XTest = Array2D.zeroCreate<float> (k+1) (k+1)
    let XTest = matMult XI X XTest

    for i = 0 to k do
         //printfn "XTest Matrix: %10.5f %10.5f " XTest.[i,0] XTest.[i,1]
         printfn "XTest Matrix: %10.5f %10.5f %10.5f %10.5f %10.5f " XTest.[i,0] XTest.[i,1] XTest.[i,2] XTest.[i,3] XTest.[i,4]

    //let a0 = Array2D.zeroCreate<float> (k+1) ni
    //let a1 = Array2D.zeroCreate<float> ni 1
    let b = matMult2 (matMult2 (matMult2 XI xp) ZI) y

    for i = 0 to k do
         printfn "Beta Hat: %i %10.5f" i b.[i,0]

    //---------------------------------------
    // Calc Std Errors
    //---------------------------------------

    let yp = transpose2 y

    let S1 = matMult2 (matMult2 yp ZI) y
    let S2 = matMult2 (matMult2 (matMult2 (transpose2 b) xp) ZI) y

    printfn "S1: %10.2f" S1.[0,0]
    printfn "S2: %10.2f" S2.[0,0]

    let sse = S1.[0,0] - S2.[0,0]
    printfn "SSE: %10.2f" sse

    let sigma2 = sse / float (ni - (k+1))
    printfn "SE: %10.2f" sigma2

    //--------------------------------------------------
    //Calc Variance/Covariance Matrix of B Hat
    //--------------------------------------------------
    //---------------------------------------
    // Calc T-Stats
    //---------------------------------------

    let bVar = Array2D.zeroCreate<float> (k+1) 1
    let bSE = Array2D.zeroCreate<float> (k+1) 1
    let TStat = Array2D.zeroCreate<float> (k+1) 1

    for i = 0 to k do
        bVar.[i,0] <- sigma2 * XI.[i,i]
        bSE.[i,0] <- Math.Sqrt bVar.[i,0]
        TStat.[i,0] <- b.[i,0] / bSE.[i,0]
        printfn "B, Var, SE, TStat: %10.5f %10.5f %10.5f %10.5f" b.[i,0] bVar.[i,0] bSE.[i,0] TStat.[i,0]

--------------------------------------------------------------------------------

When we run the GLS regression here is what we get:

---------------------------------------------------------------------
Beta Hat: 0 31086.62449
Beta Hat: 1 -65.68368
Beta Hat: 2 739.92914
Beta Hat: 3 -85.11546
Beta Hat: 4 11.52357

S1: 932945003919.98
S2: 927832435762.29
SSE: 5112568157.69
SE: 11862107.09

B, Var, SE, TStat: 31086.62449 7459397.27394 2731.18972   11.38208
B, Var, SE, TStat: -65.68368     166.51318   12.90400   -5.09018
B, Var, SE, TStat:   739.92914     268.59157   16.38876   45.14856
B, Var, SE, TStat:   -85.11546     414.70682   20.36435   -4.17963
B, Var, SE, TStat:    11.52357      20.11106    4.48454    2.56962
---------------------------------------------------------------------

My original straight OLS estimates are shown here:

---------------------------------------------------------------------
B, Var, SE, TStat: 29843.97538 3520606.31417 1876.32788   15.90552
B, Var, SE, TStat: -61.85082    74.88767    8.65377   -7.14727
B, Var, SE, TStat: 732.56913     126.54212   11.24909   65.12250
B, Var, SE, TStat: -68.85197     198.90621   14.10341   -4.88194
B, Var, SE, TStat:   10.10098      10.35856    3.21847    3.13844
---------------------------------------------------------------------

It turns out that the parameter estimates (Beta Hats) are pretty much the same between the two regressions. I suppose this is to be expected since both yield unbiased estimates of the parameters. The standard errors are higher in the GLS regression, but that is OK. The standard errors calculated under plain vanilla OLS are biased so they aren't any good anyway. The standard errors under GLS are unbiased and are the appropriate ones to use for hypothesis testing.

So, there you have F# & GLS!

PS. If anyone finds a problem, please let me know!

Where is the Economy - Railroad Traffic

2009-08-13T18:30:07Z

So, In my quest to see where the economy stands, I came across some really enlightening data. Namely railroad traffic data published by the Association of American of Railroads. They report data on carloads, by type of product, on a weekly basis. I figure that looking at the weekly data can give a good idea of where the economy is in terms of the recession. So, here's the first chart:

This chart traces the weekly carload data from September 1, 2009 to the week ending August 8th 2009. The carload data are for all products, but exclude intermodal (trailers & containers) carloads. For reference, the lowest bar is for the Christmas/New Year holidays. The trend line is a moving 4-week average. You can see that traffic really fell off during the 4th quarter of 2008 and has stayed low. The last 3 or 4 weeks show some uptick in traffic, but it would be hard to conclude that things are getting rosy fast. The recent Fed statement that the economy seems to have leveled is borne out in these data.

I also plotted Coal shipments (by far the largest product category):

Coal shipments are an interesting product to look at because they are so tied to electrical generation. We see basically the same pattern as above, with again a slight uptick in the last 3-4 weeks. I should note weather plays a role in electricity generation as a warm summer will yield extra demand due to air conditioning. I haven't looked at the data, but it seems like this summer has been cooler than last.

Next I looked at chemicals - the next largest product category:

Generally the chemicals that are shipped by rail are commodity type chemicals. We see a somewhat muted pattern for this product category -- lower than last year, but pretty steady. (I should check on that first data point, it looks out of whack)

Next we look at crushed stone, sand & gravel. Perhaps an indicator of highway and other infrastructure spending. Are the stimulus dollars having an impact?

Things still look flat.

Lastly, I took a look at intermodal traffic -- truck trailers and containers. This category is large and would seem to serve as a good indicator of how imports are faring.

Here we see the biggest drop off in the 4th quarter of 2008. Things have been flat since March, with perhaps a slight uptick recently. Nothing there to really write home about.

There you have it.

I think this data source is very valuable because it is timely -- it comes out weekly on Thursdays or the following week. Because there are only a few railroads (and all of these are members of the AAR) that make up most of the freight traffic, the survey almost covers perhaps 90% of the total rail volume.