peijun

Refer to the questions here.

Q1 Age and Voc size

dta1 <- read.table("http://140.116.183.121/~sheu/lmm/Data/age_vocab.txt", h = T)

# original plot
plot(Voc_size ~ Month, dta1)

# model
b_init <- c(b1 = 2900, b2 = -10, b3 = -0.05)
summary(m0 <- nls(Voc_size ~ b1*exp(b2*exp(b3*Month)), data = dta1, start = b_init))

Formula: Voc_size ~ b1 * exp(b2 * exp(b3 * Month))

Parameters:
     Estimate Std. Error t value Pr(>|t|)    
b1  2.905e+03  9.914e+01  29.305 1.55e-12 ***
b2 -1.018e+01  1.048e+00  -9.717 4.88e-07 ***
b3 -5.814e-02  3.446e-03 -16.873 1.00e-09 ***
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Residual standard error: 45.15 on 12 degrees of freedom

Number of iterations to convergence: 6 
Achieved convergence tolerance: 2.978e-06
#plot
ggplot(dta1, aes(Month, Voc_size, color = blue)) +
  stat_smooth(method = "nls", formula = y ~ b1*exp(b2*exp(b3*x)),
              method.args = list(start = b_init), se = F, col = 2,
              size = rel(1)) +
  geom_point(color = "blue" , size = rel(1.5)) +
  scale_x_continuous(limits = c(5, 70), breaks = seq(10, 70, by = 10)) +
  scale_y_continuous(limits = c(0, 2500), breaks = seq(0, 2500, by = 500)) +
  geom_vline(xintercept = c(40, 50), linetype = "dotted") +
  labs(x = "Age(month)", y = "Volcabulary size") +
  theme_bw()
Warning: Removed 1 rows containing non-finite values (stat_smooth).
Warning: Removed 1 rows containing missing values (geom_point).

plot(m0)

#model
b1_init <- c(b1 = 2500, b2 = 19.5, b3=-0.1)
summary(m1 <- nls(Voc_size ~ b1 * b1/(1 + (((b1 - b2)/b2) * exp(b3 * Month)))
                  , data = dta1, start = b1_init))

Formula: Voc_size ~ b1 * b1/(1 + (((b1 - b2)/b2) * exp(b3 * Month)))

Parameters:
    Estimate Std. Error t value Pr(>|t|)    
b1 50.200702   0.942465  53.265 1.26e-15 ***
b2  0.387613   0.137216   2.825   0.0153 *  
b3 -0.111637   0.009518 -11.729 6.24e-08 ***
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Residual standard error: 89.03 on 12 degrees of freedom

Number of iterations to convergence: 13 
Achieved convergence tolerance: 4.104e-06
plot(m1)

Q2 Percent of recall

dta2 <- read.table("http://140.116.183.121/~sheu/lmm/Data/freeRecall.asc", h = T)

#model
b_init <- list(a = c(5, 10), b = c(1, 0.3))
summary(m0 <- nls(ncr ~ a[grp] * exp(b[grp] * sqrt(trial)), data = dta2,
                  start = b_init))

Formula: ncr ~ a[grp] * exp(b[grp] * sqrt(trial))

Parameters:
   Estimate Std. Error t value Pr(>|t|)    
a1  7.74164    0.54481  14.210 1.72e-10 ***
a2 10.30865    0.65693  15.692 3.88e-11 ***
b1  0.23217    0.02795   8.305 3.40e-07 ***
b2  0.12121    0.02612   4.640 0.000272 ***
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Residual standard error: 0.7271 on 16 degrees of freedom

Number of iterations to convergence: 7 
Achieved convergence tolerance: 3.637e-06
#plot
dta2$yhat <- fitted(m0)
ggplot(dta2, aes(trial, ncr, color = grp)) +
  geom_point(size = 2) +
  geom_line(aes(trial,yhat), linetype = 2, size = 1) +
  labs(x = "Trials on List B", y = "Mean Numbers of correct reponses") +
  scale_color_manual(values = c("light blue", "orange")) +
  scale_x_continuous(limits = c(0, 10), breaks = seq(0, 10, by = 2)) +
  scale_y_continuous(limits = c(5, 20), breaks = seq(5, 20, by = 5)) +
  theme_bw()

Q3 Accuracy of monitoring task

dta3 <- read.table("http://140.116.183.121/~sheu/lmm/Data/continuousMonitoring.asc", h = T)

plot(accuracy ~ rate, dta3)

#model
b_init <- c(a = 20, b = 80, c = 10)
summary(m0 <- nls(accuracy ~ a + ((b-a) * rate) / (c + rate), data = dta3, start = b_init))

Formula: accuracy ~ a + ((b - a) * rate)/(c + rate)

Parameters:
  Estimate Std. Error t value Pr(>|t|)    
a   21.879      6.111   3.580  0.00128 ** 
b  100.962      8.094  12.474 5.94e-13 ***
c   10.677      4.400   2.426  0.02195 *  
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Residual standard error: 7.229 on 28 degrees of freedom

Number of iterations to convergence: 5 
Achieved convergence tolerance: 1.579e-06
summary(m1 <- nls(accuracy ~ a + b * rate + c * rate^2, data = dta3, start = b_init))

Formula: accuracy ~ a + b * rate + c * rate^2

Parameters:
  Estimate Std. Error t value Pr(>|t|)    
a 31.37523    2.92410  10.730 1.99e-11 ***
b  2.99435    0.46866   6.389 6.46e-07 ***
c -0.04255    0.01143  -3.723 0.000879 ***
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Residual standard error: 7.866 on 28 degrees of freedom

Number of iterations to convergence: 1 
Achieved convergence tolerance: 1.743e-08
# plot
yhat0 <- fitted(m0)
yhat1 <- fitted(m1)

ggplot(dta3, aes(rate, accuracy)) +
  geom_point(size = rel(2), color = "blue") +
  geom_line(aes(rate, yhat0), linetype = 2, size = .5) +
  labs(x = "Rate of stimulus change (ms/100)", y = "Monitoring accuracy (%)") +
  scale_x_continuous(limits = c(0, 45), breaks = seq(0, 45, by = 5)) +
  scale_y_continuous(limits = c(20, 100), breaks = seq(20, 100, by = 20)) +
  theme_bw()

ggplot(dta3, aes(rate, accuracy)) +
  geom_point(size = rel(2), color = "blue") +
  geom_line(aes(rate, yhat1), linetype = 2, size = .5) +
  labs(x = "Rate of stimulus change (ms/100)", y = "Monitoring accuracy (%)") +
  scale_x_continuous(limits = c(0, 45), breaks = seq(0, 45, by = 5)) +
  scale_y_continuous(limits = c(20, 100), breaks = seq(20, 100, by = 20)) +
  theme_bw()

plot(m1)

Q4 Year and population growth

dta4 <- us_national_population

ggplot(dta4, aes(year, population))+
  geom_point(pch = 1, size = 2)+
  stat_smooth(method = "gam", formula = y ~ s(x), se = F)+
  labs(x = "Year", y = "Population")+
  theme_bw()

b_init = c(b1 = 300000000 , b2 = 0, b3 = 300000000/200)
#summary(m0 <- nls(population ~ b1 / (1 + exp(b2 - b3*year)), data = dta4, start = b_init))

#Error in nlsModel(formula, mf, start, wts) : 
#  singular gradient matrix at initial parameter estimates