[Rsiena-commits] r151 - in pkg: RSiena RSiena/R RSiena/man RSienaTest RSienaTest/R RSienaTest/man
noreply at r-forge.r-project.org
noreply at r-forge.r-project.org
Fri May 27 18:17:01 CEST 2011
Author: ripleyrm
Date: 2011-05-27 18:17:01 +0200 (Fri, 27 May 2011)
New Revision: 151
Modified:
pkg/RSiena/DESCRIPTION
pkg/RSiena/R/effects.r
pkg/RSiena/changeLog
pkg/RSiena/man/RSiena-package.Rd
pkg/RSienaTest/DESCRIPTION
pkg/RSienaTest/R/effects.r
pkg/RSienaTest/R/observationErrors.r
pkg/RSienaTest/changeLog
pkg/RSienaTest/man/RSiena-package.Rd
Log:
Minor edits, traps in observation errors in RSienaTest missing from previous update.
Modified: pkg/RSiena/DESCRIPTION
===================================================================
--- pkg/RSiena/DESCRIPTION 2011-05-27 15:20:17 UTC (rev 150)
+++ pkg/RSiena/DESCRIPTION 2011-05-27 16:17:01 UTC (rev 151)
@@ -1,7 +1,7 @@
Package: RSiena
Type: Package
Title: Siena - Simulation Investigation for Empirical Network Analysis
-Version: 1.0.12.150
+Version: 1.0.12.151
Date: 2011-05-27
Author: Various
Depends: R (>= 2.10.0), xtable
Modified: pkg/RSiena/R/effects.r
===================================================================
--- pkg/RSiena/R/effects.r 2011-05-27 15:20:17 UTC (rev 150)
+++ pkg/RSiena/R/effects.r 2011-05-27 16:17:01 UTC (rev 151)
@@ -905,12 +905,12 @@
}
n <- length(xx$depvars)
types <- sapply(xx$depvars, function(x)attr(x, 'type'))
- sparses <- sapply(xx$depvars, function(x)attr(x, 'sparse'))
+ #sparses <- sapply(xx$depvars, function(x)attr(x, 'sparse'))
nOneModes <- sum(types == 'oneMode')
- nBehaviors <- sum(types == 'behavior')
+ #nBehaviors <- sum(types == 'behavior')
nBipartites <- sum(types =='bipartite')
effects <- vector('list',n)
- nodeSetNames <- sapply(xx$nodeSets, function(x)attr(x, 'nodeSetName'))
+ #nodeSetNames <- sapply(xx$nodeSets, function(x)attr(x, 'nodeSetName'))
names(effects) <- names(xx$depvars)
for (i in 1:n)
{
Modified: pkg/RSiena/changeLog
===================================================================
--- pkg/RSiena/changeLog 2011-05-27 15:20:17 UTC (rev 150)
+++ pkg/RSiena/changeLog 2011-05-27 16:17:01 UTC (rev 151)
@@ -1,3 +1,10 @@
+2011-05-27 R-forge revision 151
+
+ * R/observationErrors.r: added various traps for vectors of length
+ 1, removed spaces from function comments (RSienaTest was missed in
+ revision 150).
+ * R/effects.r: remove final few unused variables.
+
2011-05-27 R-forge revision 150
* cleanup.win, cleanup, src/Makevars, src/Makevars.win: changes
Modified: pkg/RSiena/man/RSiena-package.Rd
===================================================================
--- pkg/RSiena/man/RSiena-package.Rd 2011-05-27 15:20:17 UTC (rev 150)
+++ pkg/RSiena/man/RSiena-package.Rd 2011-05-27 16:17:01 UTC (rev 151)
@@ -30,7 +30,7 @@
\tabular{ll}{
Package: \tab RSiena\cr
Type: \tab Package\cr
-Version: \tab 1.0.12.150\cr
+Version: \tab 1.0.12.151\cr
Date: \tab 2011-05-27\cr
License: \tab GPL-2 \cr
LazyLoad: \tab yes\cr
Modified: pkg/RSienaTest/DESCRIPTION
===================================================================
--- pkg/RSienaTest/DESCRIPTION 2011-05-27 15:20:17 UTC (rev 150)
+++ pkg/RSienaTest/DESCRIPTION 2011-05-27 16:17:01 UTC (rev 151)
@@ -1,7 +1,7 @@
Package: RSienaTest
Type: Package
Title: Siena - Simulation Investigation for Empirical Network Analysis
-Version: 1.0.12.150
+Version: 1.0.12.151
Date: 2011-05-27
Author: Various
Depends: R (>= 2.10.0), xtable
Modified: pkg/RSienaTest/R/effects.r
===================================================================
--- pkg/RSienaTest/R/effects.r 2011-05-27 15:20:17 UTC (rev 150)
+++ pkg/RSienaTest/R/effects.r 2011-05-27 16:17:01 UTC (rev 151)
@@ -264,26 +264,20 @@
{
objEffects <-
rbind(objEffects,
- createEffects("covarNetNetObjective",
- otherName,
- names(xx$cCovars)[k],
- name=varname,
- groupName=groupName, group=group,
- netType=netType))
+ covarNetNetEff(otherName, names(xx$cCovars)[k],
+ attr(xx$cCovars[[k]], 'poszvar'),
+ name=varname))
}
}
- for (k in seq(along=xx$vCovars))
+ for (k in seq(along=xx$vCovars))
{
if (attr(xx$vCovars[[k]], 'nodeSet') == nodeSet)
{
objEffects <-
rbind(objEffects,
- createEffects("covarNetNetObjective",
- otherName,
- names(xx$vCovars)[k],
- name=varname,
- groupName=groupName, group=group,
- netType=netType))
+ covarNetNetEff(otherName, names(xx$vCovars)[k],
+ attr(xx$vCovars[[k]], 'poszvar'),
+ name=varname))
}
}
for (k in seq(along=xx$depvars))
@@ -293,12 +287,9 @@
{
objEffects <-
rbind(objEffects,
- createEffects("covarNetNetObjective",
- otherName,
- names(xx$depvars)[k],
- name=varname,
- groupName=groupName, group=group,
- netType=netType))
+ covarNetNetEff(otherName, names(xx$depvars)[k],
+ poszvar=TRUE,
+ name=varname))
}
}
@@ -862,13 +853,13 @@
}
##@covarNetNetEff internal getEffects
covarNetNetEff<- function(othernetname,
- covarname, poszvar, moreThan2, name)
+ covarname, poszvar, name)
{
objEffects <- createEffects("covarNetNetObjective", othernetname,
covarname, name=name,
groupName=groupName, group=group,
netType=netType)
- if (!(poszvar && (!moreThan2)))
+ if (!poszvar)
{
objEffects <- objEffects[objEffects$shortName != "covNetNet", ]
}
@@ -914,12 +905,12 @@
}
n <- length(xx$depvars)
types <- sapply(xx$depvars, function(x)attr(x, 'type'))
- sparses <- sapply(xx$depvars, function(x)attr(x, 'sparse'))
+ #sparses <- sapply(xx$depvars, function(x)attr(x, 'sparse'))
nOneModes <- sum(types == 'oneMode')
- nBehaviors <- sum(types == 'behavior')
+ #nBehaviors <- sum(types == 'behavior')
nBipartites <- sum(types =='bipartite')
effects <- vector('list',n)
- nodeSetNames <- sapply(xx$nodeSets, function(x)attr(x, 'nodeSetName'))
+ #nodeSetNames <- sapply(xx$nodeSets, function(x)attr(x, 'nodeSetName'))
names(effects) <- names(xx$depvars)
for (i in 1:n)
{
Modified: pkg/RSienaTest/R/observationErrors.r
===================================================================
--- pkg/RSienaTest/R/observationErrors.r 2011-05-27 15:20:17 UTC (rev 150)
+++ pkg/RSienaTest/R/observationErrors.r 2011-05-27 16:17:01 UTC (rev 151)
@@ -17,9 +17,9 @@
##@deviance.observationErrors siena08: deviance under normality assumptions
deviance.observationErrors <- function(mu, sig, x, se)
{
-## deviance.observationErrors: given vector observations x, assuming model
-## x ~ normal(expected = mu, variance = sig^2 + se^2)
-## required is length(x) = length(se)
+ ## deviance.observationErrors: given vector observations x, assuming model
+ ## x ~ normal(expected = mu, variance = sig^2 + se^2)
+ ## required is length(x) = length(se)
sig2 <- sig^2
sum((x - mu)^2/(sig2 + se^2)) + sum(log(sig2 + se^2))
}
@@ -27,36 +27,46 @@
##@profdev.mu siena08 (profile deviance =) minus twice profile loglik for mu
profdev.mu <- function(mu,x,se)
{
-## deviance.observationErrors maximized over sigma
+ ## deviance.observationErrors maximized over sigma
ra <- max(x) - min(x) # easy upper bound for sigma
optimize(function(sig)deviance.observationErrors(mu, sig, x, se),
- c(0, ra))
+ c(0, ra))
}
##@profdev.sig siena08 (profile deviance =) minus twice prof. loglik for sigma
profdev.sig <- function(sig,x,se)
{
-## deviance.observationErrors maximized over mu
+ ## deviance.observationErrors maximized over mu
optimize(function(mu)deviance.observationErrors(mu, sig, x, se),
- range(x))
+ range(x))
}
##@maxlik siena08 maximum likelihood estimator
maxlik <- function(x,se)
{
- ## deviance.observationErrors minimized over mu and sigma
- ## Minimize profile deviance for mu:
- opmu <- optimize(function(mu)profdev.mu(mu, x, se)$objective,
- range(x))
- ## MLE for mu:
- mu <- opmu$minimum
- ## minimized deviance:
- dev <- opmu$objective
- ## Location for minimum of deviance for this value of mu:
- sig <- profdev.mu(mu, x ,se)$minimum
- ## Standard error of MLE(mu):
- se.mu <- sqrt(1 / sum(1 / (se^2 + sig^2)))
+ if (length(x) > 1)
+ {
+ ## deviance.observationErrors minimized over mu and sigma
+ ## Minimize profile deviance for mu:
+ opmu <- optimize(function(mu)profdev.mu(mu, x, se)$objective,
+ range(x))
+ ## MLE for mu:
+ mu <- opmu$minimum
+ ## minimized deviance:
+ dev <- opmu$objective
+ ## Location for minimum of deviance for this value of mu:
+ sig <- profdev.mu(mu, x , se)$minimum
+ ## Standard error of MLE(mu):
+ se.mu <- sqrt(1 / sum(1 / (se^2 + sig^2)))
+ }
+ else
+ {
+ mu <- x
+ se.mu <- NA
+ sig <- NA
+ dev <- NA
+ }
return(list(mu = mu, se.mu = se.mu, sigma = sig, deviance = dev))
}
@@ -88,7 +98,7 @@
{
x2 <- x2 + ra
}
- if (f(x1)*f(x2) < 0)
+ if (f(x1) * f(x2) < 0)
{
break
}
@@ -104,7 +114,7 @@
}
}
-##@ confint.mu siena08 confidence interval for mu
+##@confint.mu siena08 confidence interval for mu
confint.mu <- function(x, se, alpha=0.05)
{
## confidence interval for mu in model
@@ -118,28 +128,36 @@
mindev <- maxli$deviance # minimized deviance
mlemu <- maxli$mu # MLE for mu
chidev <- qchisq(1 - alpha, 1) # critical value chi-squared distribution
- ## The end points of the confidence interval are the values
- ## where the profile deviance for mu is equal to mindev + chidev.
- ## Now first the solution for the left side of the confidence interval:
- tmp1 <- unisroot(function(mu)
- {
- tmp <- profdev.mu(mu, x, se)
- tmp$objective - mindev - chidev
- },
- c(mi, mlemu))
- mu1 <- tmp1$root
- ## and then the solution for the right side of the confidence interval:
- tmp1 <- unisroot(function(mu)
- {
- tmp <- profdev.mu(mu, x, se)
- tmp$objective - mindev - chidev
- },
- c(mlemu, ma), left=FALSE)
- mu2 <- tmp1$root
+ if (length(x) > 1)
+ {
+ ## The end points of the confidence interval are the values
+ ## where the profile deviance for mu is equal to mindev + chidev.
+ ## Now first the solution for the left side of the confidence interval:
+ tmp1 <- unisroot(function(mu)
+ {
+ tmp <- profdev.mu(mu, x, se)
+ tmp$objective - mindev - chidev
+ },
+ c(mi, mlemu))
+ mu1 <- tmp1$root
+ ## and then the solution for the right side of the confidence interval:
+ tmp1 <- unisroot(function(mu)
+ {
+ tmp <- profdev.mu(mu, x, se)
+ tmp$objective - mindev - chidev
+ },
+ c(mlemu, ma), left=FALSE)
+ mu2 <- tmp1$root
+ }
+ else
+ {
+ mu1 <- NA
+ mu2 <- NA
+ }
c(mu1, mu2, 1 - alpha)
}
-##@ confint.sig siena08 confidence interval for sigma
+##@confint.sig siena08 confidence interval for sigma
confint.sig <- function(x, se, alpha=0.05)
{
## confidence interval for sigma in model
@@ -157,20 +175,28 @@
## unless the profile deviance in 0 is smaller than this.
## Now first the solution for the left side of the confidence interval;
## this may be 0.
- if (profdev.sig(0, x, se)$objective <= mindev + chidev)
+ if (length(x) > 1)
{
- sig1 <- 0
+ if (profdev.sig(0, x, se)$objective <= mindev + chidev)
+ {
+ sig1 <- 0
+ }
+ else
+ {
+ sig1 <- uniroot(function(sig)
+ profdev.sig(sig, x, se)$objective - mindev - chidev,
+ c(0, mlesig))$root
+ }
+ ## and then the solution for the right side of the confidence interval:
+ sig2 <- unisroot(function(sig)
+ profdev.sig(sig, x, se)$objective - mindev - chidev,
+ c(mlesig, ra), left=FALSE)$root
}
else
{
- sig1 <- uniroot(function(sig)
- profdev.sig(sig, x, se)$objective - mindev - chidev,
- c(0, mlesig))$root
+ sig1 <- NA
+ sig2 <- NA
}
- ## and then the solution for the right side of the confidence interval:
- sig2 <- unisroot(function(sig)
- profdev.sig(sig,x,se)$objective - mindev - chidev,
- c(mlesig, ra), left=FALSE)$root
return(c(sig1, sig2, 1 - alpha))
}
Modified: pkg/RSienaTest/changeLog
===================================================================
--- pkg/RSienaTest/changeLog 2011-05-27 15:20:17 UTC (rev 150)
+++ pkg/RSienaTest/changeLog 2011-05-27 16:17:01 UTC (rev 151)
@@ -1,3 +1,10 @@
+2011-05-27 R-forge revision 151
+
+ * R/observationErrors.r: added various traps for vectors of length
+ 1, removed spaces from function comments (RSienaTest was missed in
+ revision 150).
+ * R/effects.r: remove final few unused variables.
+
2011-05-27 R-forge revision 150
* cleanup.win, cleanup, src/Makevars, src/Makevars.win: changes
Modified: pkg/RSienaTest/man/RSiena-package.Rd
===================================================================
--- pkg/RSienaTest/man/RSiena-package.Rd 2011-05-27 15:20:17 UTC (rev 150)
+++ pkg/RSienaTest/man/RSiena-package.Rd 2011-05-27 16:17:01 UTC (rev 151)
@@ -30,7 +30,7 @@
\tabular{ll}{
Package: \tab RSiena\cr
Type: \tab Package\cr
-Version: \tab 1.0.12.150\cr
+Version: \tab 1.0.12.151\cr
Date: \tab 2011-05-27\cr
License: \tab GPL-2 \cr
LazyLoad: \tab yes\cr
More information about the Rsiena-commits
mailing list