How do I create a link function for a power law distribution to use in glmfit?

Hi,

I am trying to use "glmfit" to generate a regression model through some data points and compare it to a simple least squared linear regression. However, I am under the assumption that my data follows a poisson distribution and obeys a power law. How do I go about creating a customized link function which can contain this information?

So far I have:

[logCoef,dev,stats]=glmfit(x, y,'poisson');

logFit=glmval(logCoef,x , 'XXX');

In the 'XXX' I want to be able to use a power law function but I have no idea how to set up this link function.

Thank you for any help and sorry if this is too vague, I am abit new to this!

Becky

Answer 1

[the cyclist]

There is documentation, and an example, of how to use a custom link in the documentation page for glmfit. Have you looked at that?

You need to be able to define the link function, its derivative, and its inverse.

Here are those functions from the example:

link = @(mu) log(mu ./ (1-mu));
derlink = @(mu) 1 ./ (mu .* (1-mu));
invlink = @(resp) 1 ./ (1 + exp(-resp));
F = {link, derlink, invlink};

and the calling syntax is

b = glmfit(X,y,'binomial','link',F)

[Rebecca]

Hi,

Thanks for this. Yes I had seen it, and I must admit I am very confused by it - I am new to writing functions!

What are mu and resp? When I put the above code into Matlab it runs, so they must be already defined by matlab?! Would I still use mu and resp when writing it for a power law?

Thanks again,

Becky

[the cyclist]

In what I copied from the documentation, link, derlink, and invlink are all "handles" to anonymous functions. This is a method of creating a user-defined function that can be used as an argument to some other function.

A simple squaring function would be written:

sqr = @(x) x.^2;

If you type that into the command window, and then type

sqr(7)

you will get 49 as an answer. You see in this case that x is a variable. It is NOT defined ahead of time. It is just a function argument.

Back to the glmfit:

link = @(mu) log(mu ./ (1-mu));

is a function definition; in this case for a logit. derlink and "invlink" are its derivative and inverse, respectively. You need to define your own link function that is a power law, and calculate the other two functions from it.

[The reason mu and response are used as variable names is just a convention referring to the mean response of the model, I presume. You can use those variable names, or not.]