Created
July 13, 2018 15:04
-
-
Save khakieconomics/2272cd7f0d61950852622198b26a2d02 to your computer and use it in GitHub Desktop.
Produce a matrix of b-splines for use internally within a Stan program
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
functions { | |
// B function | |
vector B(vector x, vector t, int i, int j_p1); | |
vector B(vector x, vector t, int i, int k) { | |
vector[rows(x)] out; | |
vector[rows(x)] a_i1j1; | |
vector[rows(x)] a_ij1; | |
if(k==1) { | |
for(n in 1:rows(x)) { | |
if((t[i] <= x[n] && x[n] < t[i+1] )){ | |
out[n] = 1.0; | |
} else { | |
out[n] = 0.0; | |
} | |
} | |
} else { | |
if(t[i+k-1] == t[i]) { | |
a_ij1 = rep_vector(0.0, rows(x)); | |
} else { | |
a_ij1 = (x - t[i]) ./ (t[i+k-1] - t[i]); | |
} | |
if(t[i+k] == t[i+1]) { | |
a_i1j1 = rep_vector(0.0, rows(x)); | |
} else { | |
a_i1j1 = (t[i+k] - x) ./ (t[i+k] - t[i+1]); | |
} | |
out = a_ij1 .* B(x, t, i, k-1) + a_i1j1 .* B(x, t, i+1, k-1); | |
} | |
return(out); | |
} | |
// Create b-spline matrix | |
matrix b_spline(vector x, vector interior_knots, int order) { | |
real xmin; | |
real xmax; | |
vector[rows(interior_knots) + 2*order] t; | |
matrix[rows(x), rows(interior_knots) + order] out; | |
// fill out augmented knots vector | |
xmin = min(x); | |
xmax = max(x); | |
t= append_row(append_row(rep_vector(xmin, order), | |
interior_knots), | |
rep_vector(xmax, order)); | |
for(p in 1:(rows(interior_knots) + order)) { | |
out[1:rows(x),p] = B(x, t, p, order-1); | |
} | |
for(i in 1:rows(out)) { | |
if(x[i]==xmax) { | |
out[i,cols(out)] = 1; | |
} | |
} | |
return(out[1:rows(out),2:cols(out)]); | |
} | |
} | |
data { | |
} | |
transformed data { | |
} | |
parameters { | |
} | |
transformed parameters { | |
} | |
model { | |
} |
Sign up for free
to join this conversation on GitHub.
Already have an account?
Sign in to comment