finterp: Formula Interpreter
In gamlss.nl: Fitting non linear parametric GAMLSS models
Description
Usage
Arguments
Value
Author(s)
References
Examples
Description
This function is taken from Jim Lindsey's R package rmutil.
What follows is taken from the help file of rmutil.
Note that not all the functionalities of finterp are implemented in nlgamlss.
finterp translates a model formula into a function of the
unknown parameters or of a vector of them. Such language formulae
can either be in Wilkinson and Rogers notation or be expressions
containing both known (existing) covariates and unknown (not
existing) parameters. In the latter, factor variables cannot be
used and parameters must be scalars.
The covariates in the formula are sought in the environment or in
the data object provided. If the data object has class, 'repeated'
or 'response', then the key words, 'times' will use the response
times from the data object as a covariate, 'individuals' will use
the index for individuals as a factor covariate, and 'nesting' the
index for nesting as a factor covariate. The latter two only work
for W&R notation.
Note that, in parameter displays, formulae in Wilkinson and Rogers
notation use variable names whereas those with unknowns use the
names of these parameters, as given in the formulae, and that the
meaning of operators (*, /, :, etc.) is different in the two
cases.
The function fmobj inspects a formula and returns a list containing the
objects referred to, with indicators as to which are unknown
parameters, covariates, factor variables, and functions.
Usage
1
2
3
4
5
6
7
8
9
10
11
## Default S3 method:
finterp(.z, .envir = parent.frame(), .formula = FALSE,
.vector = TRUE, .args = NULL, .start = 1,
.name = NULL, .expand = TRUE, .intercept = TRUE,
.old = NULL, .response = FALSE, ...)
## S3 method for class 'data.frame'
finterp(.z, .envir = NULL, .formula = FALSE, .vector = TRUE,
.args = NULL, .start = 1, .name = NULL, .expand = NULL, .intercept = TRUE,
.old = NULL, ...)
finterp(.z, ...)
fmobj(z, envir = parent.frame())
Arguments
.z
A model formula beginning with \~, either in Wilkinson and
Rogers notation or containing unknown parameters. If it
contains unknown parameters, it can have several lines so
that, for example, local variables can be assigned temporary
values. In this case, enclose the formula in curly brackets
.envir
The environment in which the formula is to be interpreted or
a data object of class, 'repeated', 'tccov', or 'tvcov'.
.formula
If TRUE and the formula is in Wilkinson and Rogers notation,
just returns the formula.
.vector
If FALSE and the formula contains unknown parameters, the
function returned has them as separate arguments. If TRUE, it
has one argument, the unknowns as a vector, unless certain
parameter names are specified in '.args'. Always TRUE if
'.envir' is a data object.
.args
If '.vector' is TRUE, names of parameters that are to be
function arguments and not included in the vector.
.start
The starting index value of the parameter vector in the
function returned when '.vector' is TRUE.
.name
Character string giving the name of the data object specified
by '.envir'. Ignored unless the latter is such an object and
only necessary when 'finterp' is called within other
functions.
.expand
If TRUE, expand functions with only time-constant covariates
to return one value per observation instead of one value per
individual. Ignored unless '.envir' is an object of class,
'repeated'.
.intercept
If W&R notation is supplied and '.intercept=F', a model
function without intercept is returned.
.old
The name of an existing object of class 'formulafn' which has
common parameters with the one being created, or a list of
such objects. Only used if '.vector'=TRUE. The value of
'.start' should ensure that there is no conflict in indexing
the vector.
.response
If TRUE, any response variable can be used in the function.
If FALSE, checks are made that the response is not also used
as a covariate.
z
A model formula beginning with ~, either in Wilkinson and
Rogers notation or containing unknown parameters.
envir
The environment in which the formula is to be interpreted.
...
for extra arguments
Value
A function, of class formulafn, of the unknown parameters or of
a vector of them is returned. Its attributes give the formula
supplied, the model function produced, the covariate names, the
parameter names, and the range of values of the index of the
parameter vector. If 'formula' is TRUE and a Wilkinson and Rogers
formula was supplied, it is simply returned instead of creating a
function.
For fmobj a list, of class 'fmobj', containing a character vector
('objects') with the names of the objects used in a formula, and
logical vectors indicating which are unknown parameters
('parameters'), covariates ('covariates'), factor variables
('factors'), and functions ('functions') is returned.
Author(s)
J.K. Lindsey
References
http://popgen.unimaas.nl/~jlindsey/index.html: Jim Lindsey web page
Examples
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
# From Jim Lindsey
x1 <- rpois(20,2)
x2 <- rnorm(20)
#
# Wilkinson and Rogers formula with three parameters
fn1 <- finterp(~x1+x2)
fn1
fn1(rep(2,3))
# the same formula with unknowns
fn2 <- finterp(~b0+b1*x1+b2*x2)
fn2
fn2(rep(2,3))
#
# nonlinear formulae with unknowns
# log link
fn2a <- finterp(~exp(b0+b1*x1+b2*x2))
fn2a
fn2a(rep(0.2,3))
# parameters common to two functions
fn2b <- finterp(~c0+c1*exp(b0+b1*x1+b2*x2), .old=fn2a, .start=4)
fn2b
# function returned also depends on values of another function
fn2c <- finterp(~fn2+c1*exp(b0+b1*x1+b2*x2), .old=fn2a,
.start=4, .args="fn2")
fn2c
args(fn2c)
fn2c(rep(0.2,4),fn2(rep(2,3)))
#
# compartment model
times <- 1:20
# exp() parameters to ensure that they are positive
fn3 <- finterp(~exp(absorption-volume)/(exp(absorption)-
exp(elimination))*(exp(-exp(elimination)*times)-
exp(-exp(absorption)*times)))
fn3
fn3(log(c(0.3,3,0.2)))
# a more efficient way
# (note that parameters do not appear in the same order)
form <- ~{
ka <- exp(absorption)
ke <- exp(elimination)
ka*exp(-volume)/(ka-ke)*(exp(-ke*times)-exp(-ka*times))}
fn3a <- finterp(form)
fn3a(log(c(0.3,0.2,3)))
#
# Poisson density
y <- rpois(20,5)
fn4 <- finterp(~mu^y*exp(-mu)/gamma(y+1))
fn4
fn4(5)
dpois(y,5)
#
# Poisson likelihood
# mean parameter
fn5 <- finterp(~-y*log(mu)+mu+lgamma(y+1),.vector=FALSE)
fn5
likefn1 <- function(p) sum(fn5(mu=p))
nlm(likefn1,p=1)
mean(y)
# canonical parameter
fn5a <- finterp(~-y*theta+exp(theta)+lgamma(y+1),.vector=FALSE)
fn5a
likefn1a <- function(p) sum(fn5a(theta=p))
nlm(likefn1a,p=1)
#
# likelihood for Poisson log linear regression
y <- rpois(20,fn2a(c(0.2,1,0.4)))
nlm(likefn1,p=1)
mean(y)
likefn2 <- function(p) sum(fn5(mu=fn2a(p)))
nlm(likefn2,p=c(1,0,0))
# or
likefn2a <- function(p) sum(fn5a(theta=fn2(p)))
nlm(likefn2a,p=c(1,0,0))
#
# likelihood for Poisson nonlinear regression
y <- rpois(20,fn3(log(c(3,0.3,0.2))))
nlm(likefn1,p=1)
mean(y)
likefn3 <- function(p) sum(fn5(mu=fn3(p)))
nlm(likefn3,p=log(c(1,0.4,0.1)))
#
# envir as data objects
# y <- matrix(rnorm(20),ncol=5)
#y[3,3] <- y[2,2] <- NA
#x1 <- 1:4
#x2 <- c("a","b","c","d")
#resp <- restovec(y)
#xx <- tcctomat(x1)
#xx2 <- tcctomat(data.frame(x1,x2))
#z1 <- matrix(rnorm(20),ncol=5)
#z2 <- matrix(rnorm(20),ncol=5)
#z3 <- matrix(rnorm(20),ncol=5)
#zz <- tvctomat(z1)
#zz <- tvctomat(z2,old=zz)
#reps <- rmna(resp, ccov=xx, tvcov=zz)
#reps2 <- rmna(resp, ccov=xx2, tvcov=zz)
#rm(y, x1, x2 , z1, z2)
#
# repeated objects
#
# time-constant covariates
# Wilkinson and Rogers notation
#form1 <- ~x1
#print(fn1 <- finterp(form1, .envir=reps))
#fn1(2:3)
#print(fn1a <- finterp(form1, .envir=xx))
#fn1a(2:3)
#form1b <- ~x1+x2
#print(fn1b <- finterp(form1b, .envir=reps2))
#fn1b(2:6)
#print(fn1c <- finterp(form1b, .envir=xx2))
#fn1c(2:6)
# with unknown parameters
#form2 <- ~a+b*x1
#print(fn2 <- finterp(form2, .envir=reps))
#fn2(2:3)
#print(fn2a <- finterp(form2, .envir=xx))
#fn2a(2:3)
#
# time-varying covariates
# Wilkinson and Rogers notation
#form3 <- ~z1+z2
#print(fn3 <- finterp(form3, .envir=reps))
#fn3(2:4)
#print(fn3a <- finterp(form3, .envir=zz))
#fn3a(2:4)
# with unknown parameters
#form4 <- ~a+b*z1+c*z2
#print(fn4 <- finterp(form4, .envir=reps))
#fn4(2:4)
#print(fn4a <- finterp(form4, .envir=zz))
#fn4a(2:4)
#
# note: lengths of x1 and z2 differ
# Wilkinson and Rogers notation
#form5 <- ~x1+z2
#print(fn5 <- finterp(form5, .envir=reps))
#fn5(2:4)
# with unknown parameters
#form6 <- ~a+b*x1+c*z2
#print(fn6 <- finterp(form6, .envir=reps))
#fn6(2:4)
#
# with times
# Wilkinson and Rogers notation
#form7 <- ~x1+z2+times
#print(fn7 <- finterp(form7, .envir=reps))
#fn7(2:5)
#form7a <- ~x1+x2+z2+times
#print(fn7a <- finterp(form7a, .envir=reps2))
#fn7a(2:8)
# with unknown parameters
#form8 <- ~a+b*x1+c*z2+e*times
#print(fn8 <- finterp(form8, .envir=reps))
#fn8(2:5)
#
# with a variable not in the data object
#form9 <- ~a+b*z1+c*z2+e*z3
#print(fn9 <- finterp(form9, .envir=reps))
#fn9(2:5)
# z3 assumed to be an unknown parameter:
#fn9(2:6)
#
# multiline formula
#form10 <- ~{
# tmp <- exp(b)
# a+tmp*z1+c*z2+d*times}
#print(fn10 <- finterp(form10, .envir=reps))
#fn10(2:5)
# for fmobj
x1 <- rpois(20,2)
x2 <- rnorm(20)
x3 <- gl(2,10)
#
# W&R formula
fmobj(~x1+x2+x3)
#
# formula with unknowns
fmobj(~b0+b1*x1+b2*x2)
#
# nonlinear formulae with unknowns
# log link
fmobj(~exp(b0+b1*x1+b2*x2))
Example output
Loading required package: gamlss
Loading required package: splines
Loading required package: gamlss.data
Attaching package: 'gamlss.data'
The following object is masked from 'package:datasets':
sleep
Loading required package: gamlss.dist
Loading required package: MASS
Loading required package: nlme
Loading required package: parallel
********** GAMLSS Version 5.1-3 **********
For more on GAMLSS look at http://www.gamlss.org/
Type gamlssNews() to see new features/changes/bug fixes.
Loading required package: survival
function (.p)
as.vector(.dm %*% .p[attr(.fna, "range")[1]:attr(.fna, "range")[2]])
<environment: 0x754f468>
attr(,"formula")
~x1 + x2
attr(,"model")
[1] "(Intercept)" "x1" "x2"
attr(,"covariates")
[1] "x1" "x2"
attr(,"parameters")
[1] "p[1]" "p[2]" "p[3]"
attr(,"range")
[1] 1 3
attr(,"class")
[1] "formulafn"
[1] 5.7096027 6.7198704 6.8691010 5.3750771 6.2252719 3.8883406
[7] 9.5302160 6.0129032 7.0595477 5.9442449 6.4538147 11.8810066
[13] 0.2926266 4.1140597 1.9826807 5.3226768 5.0047638 3.7188035
[19] 6.9162759 6.1584533
function (.p)
eval(attr(.fna, "model"))
<bytecode: 0x7b459a8>
<environment: 0x7673760>
attr(,"formula")
~b0 + b1 * x1 + b2 * x2
attr(,"model")
expression(.p[1] + .p[2] * x1 + .p[3] * x2)
attr(,"parameters")
[1] "b0" "b1" "b2"
attr(,"covariates")
[1] "x1" "x2"
attr(,"range")
[1] 1 3
attr(,"class")
[1] "formulafn"
[1] 5.7096027 6.7198704 6.8691010 5.3750771 6.2252719 3.8883406
[7] 9.5302160 6.0129032 7.0595477 5.9442449 6.4538147 11.8810066
[13] 0.2926266 4.1140597 1.9826807 5.3226768 5.0047638 3.7188035
[19] 6.9162759 6.1584533
function (.p)
eval(attr(.fna, "model"))
<bytecode: 0x7b459a8>
<environment: 0x4d712f0>
attr(,"formula")
~exp(b0 + b1 * x1 + b2 * x2)
attr(,"model")
expression(exp(.p[1] + .p[2] * x1 + .p[3] * x2))
attr(,"parameters")
[1] "b0" "b1" "b2"
attr(,"covariates")
[1] "x1" "x2"
attr(,"range")
[1] 1 3
attr(,"class")
[1] "formulafn"
[1] 1.769966 1.958124 1.987565 1.711735 1.863632 1.475260 2.593534 1.824471
[9] 2.025780 1.811988 1.906714 3.280844 1.029695 1.508938 1.219289 1.702789
[17] 1.649507 1.450459 1.996963 1.851221
function (.p)
eval(attr(.fna, "model"))
<bytecode: 0x7b459a8>
<environment: 0x3c38600>
attr(,"formula")
~c0 + c1 * exp(b0 + b1 * x1 + b2 * x2)
attr(,"model")
expression(.p[4] + .p[5] * exp(.p[1] + .p[2] * x1 + .p[3] * x2))
attr(,"parameters")
[1] "c0" "c1"
attr(,"common")
[1] "b0" "b1" "b2"
attr(,"covariates")
[1] "x1" "x2"
attr(,"range")
[1] 4 5
attr(,"class")
[1] "formulafn"
function (.p, fn2)
eval(attr(.fna, "model"))
<environment: 0x2a6d500>
attr(,"formula")
~fn2 + c1 * exp(b0 + b1 * x1 + b2 * x2)
attr(,"model")
expression(fn2 + .p[4] * exp(.p[1] + .p[2] * x1 + .p[3] * x2))
attr(,"parameters")
[1] "c1"
attr(,"common")
[1] "b0" "b1" "b2"
attr(,"covariates")
[1] "x1" "x2"
attr(,"range")
[1] 4 4
attr(,"class")
[1] "formulafn"
function (.p, fn2)
NULL
[1] 6.0635959 7.1114952 7.2666139 5.7174241 6.5979983 4.1833925
[7] 10.0489229 6.3777975 7.4647037 6.3066425 6.8351575 12.5371754
[13] 0.4985656 4.4158473 2.2265385 5.6632347 5.3346651 4.0088954
[19] 7.3156685 6.5286975
function (.p)
eval(attr(.fna, "model"))
<bytecode: 0x7b459a8>
<environment: 0x29ed470>
attr(,"formula")
~exp(absorption - volume)/(exp(absorption) - exp(elimination)) *
(exp(-exp(elimination) * times) - exp(-exp(absorption) *
times))
attr(,"model")
expression(exp(.p[1] - .p[2])/(exp(.p[1]) - exp(.p[3])) * (exp(-exp(.p[3]) *
times) - exp(-exp(.p[1]) * times)))
attr(,"parameters")
[1] "absorption" "volume" "elimination"
attr(,"covariates")
[1] "times"
attr(,"range")
[1] 1 3
attr(,"class")
[1] "formulafn"
[1] 0.07791253 0.12150841 0.14224198 0.14813475 0.14474928 0.13589532
[7] 0.12414054 0.11117856 0.09809338 0.08554821 0.07391999 0.06339423
[13] 0.05403167 0.04581449 0.03867807 0.03253246 0.02727652 0.02280714
[19] 0.01902481 0.01583689
[1] 0.07791253 0.12150841 0.14224198 0.14813475 0.14474928 0.13589532
[7] 0.12414054 0.11117856 0.09809338 0.08554821 0.07391999 0.06339423
[13] 0.05403167 0.04581449 0.03867807 0.03253246 0.02727652 0.02280714
[19] 0.01902481 0.01583689
function (.p)
eval(attr(.fna, "model"))
<bytecode: 0x7b459a8>
<environment: 0x69770d8>
attr(,"formula")
~mu^y * exp(-mu)/gamma(y + 1)
attr(,"model")
expression(.p[1]^y * exp(-.p[1])/gamma(y + 1))
attr(,"parameters")
[1] "mu"
attr(,"covariates")
[1] "y"
attr(,"range")
[1] 1 1
attr(,"class")
[1] "formulafn"
[1] 0.036265577 0.084224337 0.175467370 0.140373896 0.175467370 0.033689735
[7] 0.084224337 0.175467370 0.084224337 0.146222808 0.175467370 0.033689735
[13] 0.175467370 0.006737947 0.140373896 0.175467370 0.065278039 0.146222808
[19] 0.104444863 0.036265577
[1] 0.036265577 0.084224337 0.175467370 0.140373896 0.175467370 0.033689735
[7] 0.084224337 0.175467370 0.084224337 0.146222808 0.175467370 0.033689735
[13] 0.175467370 0.006737947 0.140373896 0.175467370 0.065278039 0.146222808
[19] 0.104444863 0.036265577
function (mu)
eval(attr(.fna, "model"))
<environment: 0x69a1020>
attr(,"formula")
~-y * log(mu) + mu + lgamma(y + 1)
attr(,"model")
expression(-y * log(mu) + mu + lgamma(y + 1))
attr(,"parameters")
[1] "mu"
attr(,"covariates")
[1] "y"
attr(,"range")
[1] 1 1
attr(,"class")
[1] "formulafn"
$minimum
[1] 48.01063
$estimate
[1] 4.3
$gradient
[1] 4.62679e-10
$code
[1] 1
$iterations
[1] 11
Warning messages:
1: In log(mu) : NaNs produced
2: In nlm(likefn1, p = 1) : NA/Inf replaced by maximum positive value
3: In log(mu) : NaNs produced
4: In nlm(likefn1, p = 1) : NA/Inf replaced by maximum positive value
[1] 4.3
function (theta)
eval(attr(.fna, "model"))
<environment: 0x6ef1b30>
attr(,"formula")
~-y * theta + exp(theta) + lgamma(y + 1)
attr(,"model")
expression(-y * theta + exp(theta) + lgamma(y + 1))
attr(,"parameters")
[1] "theta"
attr(,"covariates")
[1] "y"
attr(,"range")
[1] 1 1
attr(,"class")
[1] "formulafn"
$minimum
[1] 48.01063
$estimate
[1] 1.458614
$gradient
[1] -5.845626e-06
$code
[1] 1
$iterations
[1] 5
$minimum
[1] 207.841
$estimate
[1] 13.34999
$gradient
[1] 1.173061e-06
$code
[1] 1
$iterations
[1] 12
Warning messages:
1: In log(mu) : NaNs produced
2: In nlm(likefn1, p = 1) : NA/Inf replaced by maximum positive value
3: In log(mu) : NaNs produced
4: In nlm(likefn1, p = 1) : NA/Inf replaced by maximum positive value
5: In log(mu) : NaNs produced
6: In nlm(likefn1, p = 1) : NA/Inf replaced by maximum positive value
[1] 13.35
$minimum
[1] 45.03136
$estimate
[1] 0.1633968 0.9846741 0.4336110
$gradient
[1] 1.197265e-05 2.921752e-05 3.673506e-06
$code
[1] 1
$iterations
[1] 14
Warning message:
In nlm(likefn2, p = c(1, 0, 0)) : NA/Inf replaced by maximum positive value
$minimum
[1] 45.03136
$estimate
[1] 0.1633968 0.9846741 0.4336110
$gradient
[1] 1.197265e-05 2.921041e-05 3.673506e-06
$code
[1] 1
$iterations
[1] 14
Warning message:
In nlm(likefn2a, p = c(1, 0, 0)) :
NA/Inf replaced by maximum positive value
$minimum
[1] 23.18643
$estimate
[1] 0.7499995
$gradient
[1] 3.161915e-07
$code
[1] 1
$iterations
[1] 4
Warning messages:
1: In log(mu) : NaNs produced
2: In nlm(likefn1, p = 1) : NA/Inf replaced by maximum positive value
[1] 0.75
$minimum
[1] 18.89604
$estimate
[1] 3.1175446 -0.8944346 -1.9378504
$gradient
[1] -3.190844e-07 -4.199663e-05 -3.764371e-05
$code
[1] 2
$iterations
[1] 34
$formula
[1] " \n x1 + x2 + x3 "
$objects
[1] "x1" "x2" "x3"
$functions
[1] "~" "+"
$parameters
[1] FALSE FALSE FALSE
$covariates
[1] TRUE TRUE FALSE
$factors
[1] FALSE FALSE TRUE
$local
[1] FALSE FALSE FALSE
attr(,"class")
[1] "fmobj"
$formula
[1] " \n b0 + b1 * x1 + b2 * x2 "
$objects
[1] "b0" "b1" "x1" "b2" "x2"
$functions
[1] "~" "+" "*"
$parameters
[1] TRUE TRUE FALSE TRUE FALSE
$covariates
[1] FALSE FALSE TRUE FALSE TRUE
$factors
[1] FALSE FALSE FALSE FALSE FALSE
$local
[1] FALSE FALSE FALSE FALSE FALSE
attr(,"class")
[1] "fmobj"
$formula
[1] " \n exp ( b0 + b1 * x1 + b2 * x2 ) "
$objects
[1] "b0" "b1" "x1" "b2" "x2"
$functions
[1] "~" "exp" "+" "*"
$parameters
[1] TRUE TRUE FALSE TRUE FALSE
$covariates
[1] FALSE FALSE TRUE FALSE TRUE
$factors
[1] FALSE FALSE FALSE FALSE FALSE
$local
[1] FALSE FALSE FALSE FALSE FALSE
attr(,"class")
[1] "fmobj"
gamlss.nl documentation built on May 2, 2019, 9:27 a.m.
R Package Documentation
Browse R Packages
We want your feedback!
Note that we can't provide technical support on individual packages. You should contact the package authors for that.
Description Usage Arguments Value Author(s) References Examples
Description
This function is taken from Jim Lindsey's R package rmutil.
What follows is taken from the help file of rmutil.
Note that not all the functionalities of finterp are implemented in nlgamlss.
finterp translates a model formula into a function of the
unknown parameters or of a vector of them. Such language formulae
can either be in Wilkinson and Rogers notation or be expressions
containing both known (existing) covariates and unknown (not
existing) parameters. In the latter, factor variables cannot be
used and parameters must be scalars.
The covariates in the formula are sought in the environment or in the data object provided. If the data object has class, 'repeated' or 'response', then the key words, 'times' will use the response times from the data object as a covariate, 'individuals' will use the index for individuals as a factor covariate, and 'nesting' the index for nesting as a factor covariate. The latter two only work for W&R notation.
Note that, in parameter displays, formulae in Wilkinson and Rogers notation use variable names whereas those with unknowns use the names of these parameters, as given in the formulae, and that the meaning of operators (*, /, :, etc.) is different in the two cases.
The function fmobj inspects a formula and returns a list containing the
objects referred to, with indicators as to which are unknown
parameters, covariates, factor variables, and functions.
Usage
1 2 3 4 5 6 7 8 9 10 11 | ## Default S3 method:
finterp(.z, .envir = parent.frame(), .formula = FALSE,
.vector = TRUE, .args = NULL, .start = 1,
.name = NULL, .expand = TRUE, .intercept = TRUE,
.old = NULL, .response = FALSE, ...)
## S3 method for class 'data.frame'
finterp(.z, .envir = NULL, .formula = FALSE, .vector = TRUE,
.args = NULL, .start = 1, .name = NULL, .expand = NULL, .intercept = TRUE,
.old = NULL, ...)
finterp(.z, ...)
fmobj(z, envir = parent.frame())
|
Arguments
.z |
A model formula beginning with \~, either in Wilkinson and Rogers notation or containing unknown parameters. If it contains unknown parameters, it can have several lines so that, for example, local variables can be assigned temporary values. In this case, enclose the formula in curly brackets |
.envir |
The environment in which the formula is to be interpreted or a data object of class, 'repeated', 'tccov', or 'tvcov'. |
.formula |
If TRUE and the formula is in Wilkinson and Rogers notation, just returns the formula. |
.vector |
If FALSE and the formula contains unknown parameters, the function returned has them as separate arguments. If TRUE, it has one argument, the unknowns as a vector, unless certain parameter names are specified in '.args'. Always TRUE if '.envir' is a data object. |
.args |
If '.vector' is TRUE, names of parameters that are to be function arguments and not included in the vector. |
.start |
The starting index value of the parameter vector in the function returned when '.vector' is TRUE. |
.name |
Character string giving the name of the data object specified by '.envir'. Ignored unless the latter is such an object and only necessary when 'finterp' is called within other functions. |
.expand |
If TRUE, expand functions with only time-constant covariates to return one value per observation instead of one value per individual. Ignored unless '.envir' is an object of class, 'repeated'. |
.intercept |
If W&R notation is supplied and '.intercept=F', a model function without intercept is returned. |
.old |
The name of an existing object of class 'formulafn' which has common parameters with the one being created, or a list of such objects. Only used if '.vector'=TRUE. The value of '.start' should ensure that there is no conflict in indexing the vector. |
.response |
If TRUE, any response variable can be used in the function. If FALSE, checks are made that the response is not also used as a covariate. |
z |
A model formula beginning with ~, either in Wilkinson and Rogers notation or containing unknown parameters. |
envir |
The environment in which the formula is to be interpreted. |
... |
for extra arguments |
Value
A function, of class formulafn, of the unknown parameters or of
a vector of them is returned. Its attributes give the formula
supplied, the model function produced, the covariate names, the
parameter names, and the range of values of the index of the
parameter vector. If 'formula' is TRUE and a Wilkinson and Rogers
formula was supplied, it is simply returned instead of creating a
function.
For fmobj a list, of class 'fmobj', containing a character vector
('objects') with the names of the objects used in a formula, and
logical vectors indicating which are unknown parameters
('parameters'), covariates ('covariates'), factor variables
('factors'), and functions ('functions') is returned.
Author(s)
J.K. Lindsey
References
http://popgen.unimaas.nl/~jlindsey/index.html: Jim Lindsey web page
Examples
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 90 91 92 93 94 95 96 97 98 99 100 101 102 103 104 105 106 107 108 109 110 111 112 113 114 115 116 117 118 119 120 121 122 123 124 125 126 127 128 129 130 131 132 133 134 135 136 137 138 139 140 141 142 143 144 145 146 147 148 149 150 151 152 153 154 155 156 157 158 159 160 161 162 163 164 165 166 167 168 169 170 171 172 173 174 175 176 177 178 179 180 181 182 183 184 | # From Jim Lindsey
x1 <- rpois(20,2)
x2 <- rnorm(20)
#
# Wilkinson and Rogers formula with three parameters
fn1 <- finterp(~x1+x2)
fn1
fn1(rep(2,3))
# the same formula with unknowns
fn2 <- finterp(~b0+b1*x1+b2*x2)
fn2
fn2(rep(2,3))
#
# nonlinear formulae with unknowns
# log link
fn2a <- finterp(~exp(b0+b1*x1+b2*x2))
fn2a
fn2a(rep(0.2,3))
# parameters common to two functions
fn2b <- finterp(~c0+c1*exp(b0+b1*x1+b2*x2), .old=fn2a, .start=4)
fn2b
# function returned also depends on values of another function
fn2c <- finterp(~fn2+c1*exp(b0+b1*x1+b2*x2), .old=fn2a,
.start=4, .args="fn2")
fn2c
args(fn2c)
fn2c(rep(0.2,4),fn2(rep(2,3)))
#
# compartment model
times <- 1:20
# exp() parameters to ensure that they are positive
fn3 <- finterp(~exp(absorption-volume)/(exp(absorption)-
exp(elimination))*(exp(-exp(elimination)*times)-
exp(-exp(absorption)*times)))
fn3
fn3(log(c(0.3,3,0.2)))
# a more efficient way
# (note that parameters do not appear in the same order)
form <- ~{
ka <- exp(absorption)
ke <- exp(elimination)
ka*exp(-volume)/(ka-ke)*(exp(-ke*times)-exp(-ka*times))}
fn3a <- finterp(form)
fn3a(log(c(0.3,0.2,3)))
#
# Poisson density
y <- rpois(20,5)
fn4 <- finterp(~mu^y*exp(-mu)/gamma(y+1))
fn4
fn4(5)
dpois(y,5)
#
# Poisson likelihood
# mean parameter
fn5 <- finterp(~-y*log(mu)+mu+lgamma(y+1),.vector=FALSE)
fn5
likefn1 <- function(p) sum(fn5(mu=p))
nlm(likefn1,p=1)
mean(y)
# canonical parameter
fn5a <- finterp(~-y*theta+exp(theta)+lgamma(y+1),.vector=FALSE)
fn5a
likefn1a <- function(p) sum(fn5a(theta=p))
nlm(likefn1a,p=1)
#
# likelihood for Poisson log linear regression
y <- rpois(20,fn2a(c(0.2,1,0.4)))
nlm(likefn1,p=1)
mean(y)
likefn2 <- function(p) sum(fn5(mu=fn2a(p)))
nlm(likefn2,p=c(1,0,0))
# or
likefn2a <- function(p) sum(fn5a(theta=fn2(p)))
nlm(likefn2a,p=c(1,0,0))
#
# likelihood for Poisson nonlinear regression
y <- rpois(20,fn3(log(c(3,0.3,0.2))))
nlm(likefn1,p=1)
mean(y)
likefn3 <- function(p) sum(fn5(mu=fn3(p)))
nlm(likefn3,p=log(c(1,0.4,0.1)))
#
# envir as data objects
# y <- matrix(rnorm(20),ncol=5)
#y[3,3] <- y[2,2] <- NA
#x1 <- 1:4
#x2 <- c("a","b","c","d")
#resp <- restovec(y)
#xx <- tcctomat(x1)
#xx2 <- tcctomat(data.frame(x1,x2))
#z1 <- matrix(rnorm(20),ncol=5)
#z2 <- matrix(rnorm(20),ncol=5)
#z3 <- matrix(rnorm(20),ncol=5)
#zz <- tvctomat(z1)
#zz <- tvctomat(z2,old=zz)
#reps <- rmna(resp, ccov=xx, tvcov=zz)
#reps2 <- rmna(resp, ccov=xx2, tvcov=zz)
#rm(y, x1, x2 , z1, z2)
#
# repeated objects
#
# time-constant covariates
# Wilkinson and Rogers notation
#form1 <- ~x1
#print(fn1 <- finterp(form1, .envir=reps))
#fn1(2:3)
#print(fn1a <- finterp(form1, .envir=xx))
#fn1a(2:3)
#form1b <- ~x1+x2
#print(fn1b <- finterp(form1b, .envir=reps2))
#fn1b(2:6)
#print(fn1c <- finterp(form1b, .envir=xx2))
#fn1c(2:6)
# with unknown parameters
#form2 <- ~a+b*x1
#print(fn2 <- finterp(form2, .envir=reps))
#fn2(2:3)
#print(fn2a <- finterp(form2, .envir=xx))
#fn2a(2:3)
#
# time-varying covariates
# Wilkinson and Rogers notation
#form3 <- ~z1+z2
#print(fn3 <- finterp(form3, .envir=reps))
#fn3(2:4)
#print(fn3a <- finterp(form3, .envir=zz))
#fn3a(2:4)
# with unknown parameters
#form4 <- ~a+b*z1+c*z2
#print(fn4 <- finterp(form4, .envir=reps))
#fn4(2:4)
#print(fn4a <- finterp(form4, .envir=zz))
#fn4a(2:4)
#
# note: lengths of x1 and z2 differ
# Wilkinson and Rogers notation
#form5 <- ~x1+z2
#print(fn5 <- finterp(form5, .envir=reps))
#fn5(2:4)
# with unknown parameters
#form6 <- ~a+b*x1+c*z2
#print(fn6 <- finterp(form6, .envir=reps))
#fn6(2:4)
#
# with times
# Wilkinson and Rogers notation
#form7 <- ~x1+z2+times
#print(fn7 <- finterp(form7, .envir=reps))
#fn7(2:5)
#form7a <- ~x1+x2+z2+times
#print(fn7a <- finterp(form7a, .envir=reps2))
#fn7a(2:8)
# with unknown parameters
#form8 <- ~a+b*x1+c*z2+e*times
#print(fn8 <- finterp(form8, .envir=reps))
#fn8(2:5)
#
# with a variable not in the data object
#form9 <- ~a+b*z1+c*z2+e*z3
#print(fn9 <- finterp(form9, .envir=reps))
#fn9(2:5)
# z3 assumed to be an unknown parameter:
#fn9(2:6)
#
# multiline formula
#form10 <- ~{
# tmp <- exp(b)
# a+tmp*z1+c*z2+d*times}
#print(fn10 <- finterp(form10, .envir=reps))
#fn10(2:5)
# for fmobj
x1 <- rpois(20,2)
x2 <- rnorm(20)
x3 <- gl(2,10)
#
# W&R formula
fmobj(~x1+x2+x3)
#
# formula with unknowns
fmobj(~b0+b1*x1+b2*x2)
#
# nonlinear formulae with unknowns
# log link
fmobj(~exp(b0+b1*x1+b2*x2))
|
Example output
Loading required package: gamlss
Loading required package: splines
Loading required package: gamlss.data
Attaching package: 'gamlss.data'
The following object is masked from 'package:datasets':
sleep
Loading required package: gamlss.dist
Loading required package: MASS
Loading required package: nlme
Loading required package: parallel
********** GAMLSS Version 5.1-3 **********
For more on GAMLSS look at http://www.gamlss.org/
Type gamlssNews() to see new features/changes/bug fixes.
Loading required package: survival
function (.p)
as.vector(.dm %*% .p[attr(.fna, "range")[1]:attr(.fna, "range")[2]])
<environment: 0x754f468>
attr(,"formula")
~x1 + x2
attr(,"model")
[1] "(Intercept)" "x1" "x2"
attr(,"covariates")
[1] "x1" "x2"
attr(,"parameters")
[1] "p[1]" "p[2]" "p[3]"
attr(,"range")
[1] 1 3
attr(,"class")
[1] "formulafn"
[1] 5.7096027 6.7198704 6.8691010 5.3750771 6.2252719 3.8883406
[7] 9.5302160 6.0129032 7.0595477 5.9442449 6.4538147 11.8810066
[13] 0.2926266 4.1140597 1.9826807 5.3226768 5.0047638 3.7188035
[19] 6.9162759 6.1584533
function (.p)
eval(attr(.fna, "model"))
<bytecode: 0x7b459a8>
<environment: 0x7673760>
attr(,"formula")
~b0 + b1 * x1 + b2 * x2
attr(,"model")
expression(.p[1] + .p[2] * x1 + .p[3] * x2)
attr(,"parameters")
[1] "b0" "b1" "b2"
attr(,"covariates")
[1] "x1" "x2"
attr(,"range")
[1] 1 3
attr(,"class")
[1] "formulafn"
[1] 5.7096027 6.7198704 6.8691010 5.3750771 6.2252719 3.8883406
[7] 9.5302160 6.0129032 7.0595477 5.9442449 6.4538147 11.8810066
[13] 0.2926266 4.1140597 1.9826807 5.3226768 5.0047638 3.7188035
[19] 6.9162759 6.1584533
function (.p)
eval(attr(.fna, "model"))
<bytecode: 0x7b459a8>
<environment: 0x4d712f0>
attr(,"formula")
~exp(b0 + b1 * x1 + b2 * x2)
attr(,"model")
expression(exp(.p[1] + .p[2] * x1 + .p[3] * x2))
attr(,"parameters")
[1] "b0" "b1" "b2"
attr(,"covariates")
[1] "x1" "x2"
attr(,"range")
[1] 1 3
attr(,"class")
[1] "formulafn"
[1] 1.769966 1.958124 1.987565 1.711735 1.863632 1.475260 2.593534 1.824471
[9] 2.025780 1.811988 1.906714 3.280844 1.029695 1.508938 1.219289 1.702789
[17] 1.649507 1.450459 1.996963 1.851221
function (.p)
eval(attr(.fna, "model"))
<bytecode: 0x7b459a8>
<environment: 0x3c38600>
attr(,"formula")
~c0 + c1 * exp(b0 + b1 * x1 + b2 * x2)
attr(,"model")
expression(.p[4] + .p[5] * exp(.p[1] + .p[2] * x1 + .p[3] * x2))
attr(,"parameters")
[1] "c0" "c1"
attr(,"common")
[1] "b0" "b1" "b2"
attr(,"covariates")
[1] "x1" "x2"
attr(,"range")
[1] 4 5
attr(,"class")
[1] "formulafn"
function (.p, fn2)
eval(attr(.fna, "model"))
<environment: 0x2a6d500>
attr(,"formula")
~fn2 + c1 * exp(b0 + b1 * x1 + b2 * x2)
attr(,"model")
expression(fn2 + .p[4] * exp(.p[1] + .p[2] * x1 + .p[3] * x2))
attr(,"parameters")
[1] "c1"
attr(,"common")
[1] "b0" "b1" "b2"
attr(,"covariates")
[1] "x1" "x2"
attr(,"range")
[1] 4 4
attr(,"class")
[1] "formulafn"
function (.p, fn2)
NULL
[1] 6.0635959 7.1114952 7.2666139 5.7174241 6.5979983 4.1833925
[7] 10.0489229 6.3777975 7.4647037 6.3066425 6.8351575 12.5371754
[13] 0.4985656 4.4158473 2.2265385 5.6632347 5.3346651 4.0088954
[19] 7.3156685 6.5286975
function (.p)
eval(attr(.fna, "model"))
<bytecode: 0x7b459a8>
<environment: 0x29ed470>
attr(,"formula")
~exp(absorption - volume)/(exp(absorption) - exp(elimination)) *
(exp(-exp(elimination) * times) - exp(-exp(absorption) *
times))
attr(,"model")
expression(exp(.p[1] - .p[2])/(exp(.p[1]) - exp(.p[3])) * (exp(-exp(.p[3]) *
times) - exp(-exp(.p[1]) * times)))
attr(,"parameters")
[1] "absorption" "volume" "elimination"
attr(,"covariates")
[1] "times"
attr(,"range")
[1] 1 3
attr(,"class")
[1] "formulafn"
[1] 0.07791253 0.12150841 0.14224198 0.14813475 0.14474928 0.13589532
[7] 0.12414054 0.11117856 0.09809338 0.08554821 0.07391999 0.06339423
[13] 0.05403167 0.04581449 0.03867807 0.03253246 0.02727652 0.02280714
[19] 0.01902481 0.01583689
[1] 0.07791253 0.12150841 0.14224198 0.14813475 0.14474928 0.13589532
[7] 0.12414054 0.11117856 0.09809338 0.08554821 0.07391999 0.06339423
[13] 0.05403167 0.04581449 0.03867807 0.03253246 0.02727652 0.02280714
[19] 0.01902481 0.01583689
function (.p)
eval(attr(.fna, "model"))
<bytecode: 0x7b459a8>
<environment: 0x69770d8>
attr(,"formula")
~mu^y * exp(-mu)/gamma(y + 1)
attr(,"model")
expression(.p[1]^y * exp(-.p[1])/gamma(y + 1))
attr(,"parameters")
[1] "mu"
attr(,"covariates")
[1] "y"
attr(,"range")
[1] 1 1
attr(,"class")
[1] "formulafn"
[1] 0.036265577 0.084224337 0.175467370 0.140373896 0.175467370 0.033689735
[7] 0.084224337 0.175467370 0.084224337 0.146222808 0.175467370 0.033689735
[13] 0.175467370 0.006737947 0.140373896 0.175467370 0.065278039 0.146222808
[19] 0.104444863 0.036265577
[1] 0.036265577 0.084224337 0.175467370 0.140373896 0.175467370 0.033689735
[7] 0.084224337 0.175467370 0.084224337 0.146222808 0.175467370 0.033689735
[13] 0.175467370 0.006737947 0.140373896 0.175467370 0.065278039 0.146222808
[19] 0.104444863 0.036265577
function (mu)
eval(attr(.fna, "model"))
<environment: 0x69a1020>
attr(,"formula")
~-y * log(mu) + mu + lgamma(y + 1)
attr(,"model")
expression(-y * log(mu) + mu + lgamma(y + 1))
attr(,"parameters")
[1] "mu"
attr(,"covariates")
[1] "y"
attr(,"range")
[1] 1 1
attr(,"class")
[1] "formulafn"
$minimum
[1] 48.01063
$estimate
[1] 4.3
$gradient
[1] 4.62679e-10
$code
[1] 1
$iterations
[1] 11
Warning messages:
1: In log(mu) : NaNs produced
2: In nlm(likefn1, p = 1) : NA/Inf replaced by maximum positive value
3: In log(mu) : NaNs produced
4: In nlm(likefn1, p = 1) : NA/Inf replaced by maximum positive value
[1] 4.3
function (theta)
eval(attr(.fna, "model"))
<environment: 0x6ef1b30>
attr(,"formula")
~-y * theta + exp(theta) + lgamma(y + 1)
attr(,"model")
expression(-y * theta + exp(theta) + lgamma(y + 1))
attr(,"parameters")
[1] "theta"
attr(,"covariates")
[1] "y"
attr(,"range")
[1] 1 1
attr(,"class")
[1] "formulafn"
$minimum
[1] 48.01063
$estimate
[1] 1.458614
$gradient
[1] -5.845626e-06
$code
[1] 1
$iterations
[1] 5
$minimum
[1] 207.841
$estimate
[1] 13.34999
$gradient
[1] 1.173061e-06
$code
[1] 1
$iterations
[1] 12
Warning messages:
1: In log(mu) : NaNs produced
2: In nlm(likefn1, p = 1) : NA/Inf replaced by maximum positive value
3: In log(mu) : NaNs produced
4: In nlm(likefn1, p = 1) : NA/Inf replaced by maximum positive value
5: In log(mu) : NaNs produced
6: In nlm(likefn1, p = 1) : NA/Inf replaced by maximum positive value
[1] 13.35
$minimum
[1] 45.03136
$estimate
[1] 0.1633968 0.9846741 0.4336110
$gradient
[1] 1.197265e-05 2.921752e-05 3.673506e-06
$code
[1] 1
$iterations
[1] 14
Warning message:
In nlm(likefn2, p = c(1, 0, 0)) : NA/Inf replaced by maximum positive value
$minimum
[1] 45.03136
$estimate
[1] 0.1633968 0.9846741 0.4336110
$gradient
[1] 1.197265e-05 2.921041e-05 3.673506e-06
$code
[1] 1
$iterations
[1] 14
Warning message:
In nlm(likefn2a, p = c(1, 0, 0)) :
NA/Inf replaced by maximum positive value
$minimum
[1] 23.18643
$estimate
[1] 0.7499995
$gradient
[1] 3.161915e-07
$code
[1] 1
$iterations
[1] 4
Warning messages:
1: In log(mu) : NaNs produced
2: In nlm(likefn1, p = 1) : NA/Inf replaced by maximum positive value
[1] 0.75
$minimum
[1] 18.89604
$estimate
[1] 3.1175446 -0.8944346 -1.9378504
$gradient
[1] -3.190844e-07 -4.199663e-05 -3.764371e-05
$code
[1] 2
$iterations
[1] 34
$formula
[1] " \n x1 + x2 + x3 "
$objects
[1] "x1" "x2" "x3"
$functions
[1] "~" "+"
$parameters
[1] FALSE FALSE FALSE
$covariates
[1] TRUE TRUE FALSE
$factors
[1] FALSE FALSE TRUE
$local
[1] FALSE FALSE FALSE
attr(,"class")
[1] "fmobj"
$formula
[1] " \n b0 + b1 * x1 + b2 * x2 "
$objects
[1] "b0" "b1" "x1" "b2" "x2"
$functions
[1] "~" "+" "*"
$parameters
[1] TRUE TRUE FALSE TRUE FALSE
$covariates
[1] FALSE FALSE TRUE FALSE TRUE
$factors
[1] FALSE FALSE FALSE FALSE FALSE
$local
[1] FALSE FALSE FALSE FALSE FALSE
attr(,"class")
[1] "fmobj"
$formula
[1] " \n exp ( b0 + b1 * x1 + b2 * x2 ) "
$objects
[1] "b0" "b1" "x1" "b2" "x2"
$functions
[1] "~" "exp" "+" "*"
$parameters
[1] TRUE TRUE FALSE TRUE FALSE
$covariates
[1] FALSE FALSE TRUE FALSE TRUE
$factors
[1] FALSE FALSE FALSE FALSE FALSE
$local
[1] FALSE FALSE FALSE FALSE FALSE
attr(,"class")
[1] "fmobj"
R Package Documentation
Browse R Packages
We want your feedback!
Note that we can't provide technical support on individual packages. You should contact the package authors for that.
Embedding an R snippet on your website
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.