Guest User

Untitled

a guest
Jul 15th, 2018
99
0
Never
Not a member of Pastebin yet? Sign Up, it unlocks many cool features!
text 10.90 KB | None | 0 0
  1. dTMP<-structure(list(jday = c(188L, 188L, 189L, 189L, 190L, 190L, 191L,
  2. 191L, 192L, 192L, 193L, 193L, 194L, 194L, 195L, 195L, 196L, 196L,
  3. 197L, 197L, 198L, 198L, 199L, 199L, 200L, 200L, 201L, 201L, 202L,
  4. 202L, 203L, 203L, 204L, 204L, 205L, 205L, 206L, 206L, 207L, 207L,
  5. 208L, 208L, 209L, 209L, 210L, 210L, 211L, 211L, 212L, 212L, 213L,
  6. 213L, 214L, 214L, 215L, 215L, 216L, 216L, 217L, 217L, 218L, 218L,
  7. 219L, 219L, 220L, 220L, 221L, 221L, 222L, 222L, 223L, 223L, 224L,
  8. 224L, 225L, 225L, 226L, 226L, 227L, 227L, 228L, 228L, 229L, 229L,
  9. 230L, 230L, 231L, 231L, 232L, 232L, 233L, 233L, 234L, 234L, 235L,
  10. 235L, 236L, 236L, 237L, 237L, 238L, 238L, 239L, 239L, 240L, 240L,
  11. 241L, 241L, 242L, 242L, 243L, 243L, 244L, 244L, 245L, 245L, 246L,
  12. 246L, 247L, 247L, 248L, 248L, 249L, 249L, 250L, 250L, 251L, 251L,
  13. 252L, 252L, 253L, 253L, 254L, 254L, 255L, 255L, 256L, 256L, 257L,
  14. 257L, 258L, 258L, 259L, 259L, 260L, 260L, 261L, 261L, 262L, 262L,
  15. 263L, 263L, 264L, 264L, 265L, 265L, 266L, 266L, 267L, 267L, 268L,
  16. 268L, 269L, 269L, 270L, 270L, 271L, 271L, 272L, 272L, 273L, 273L,
  17. 274L, 274L, 275L, 275L, 276L, 276L, 277L, 277L, 278L, 278L, 279L,
  18. 279L, 280L, 280L, 281L, 281L), site = c("A", "B", "A", "B", "A",
  19. "B", "A", "B", "A", "B", "A", "B", "A", "B", "A", "B", "A", "B",
  20. "A", "B", "A", "B", "A", "B", "A", "B", "A", "B", "A", "B", "A",
  21. "B", "A", "B", "A", "B", "A", "B", "A", "B", "A", "B", "A", "B",
  22. "A", "B", "A", "B", "A", "B", "A", "B", "A", "B", "A", "B", "A",
  23. "B", "A", "B", "A", "B", "A", "B", "A", "B", "A", "B", "A", "B",
  24. "A", "B", "A", "B", "A", "B", "A", "B", "A", "B", "A", "B", "A",
  25. "B", "A", "B", "A", "B", "A", "B", "A", "B", "A", "B", "A", "B",
  26. "A", "B", "A", "B", "A", "B", "A", "B", "A", "B", "A", "B", "A",
  27. "B", "A", "B", "A", "B", "A", "B", "A", "B", "A", "B", "A", "B",
  28. "A", "B", "A", "B", "A", "B", "A", "B", "A", "B", "A", "B", "A",
  29. "B", "A", "B", "A", "B", "A", "B", "A", "B", "A", "B", "A", "B",
  30. "A", "B", "A", "B", "A", "B", "A", "B", "A", "B", "A", "B", "A",
  31. "B", "A", "B", "A", "B", "A", "B", "A", "B", "A", "B", "A", "B",
  32. "A", "B", "A", "B", "A", "B", "A", "B", "A", "B", "A", "B", "A",
  33. "B"), temp.avg = c(5.99097902097902, 6.09086713286713, 7.42745833333333,
  34. 8.1678125, 6.88500694444444, 8.23406944444444, 6.12922222222222,
  35. 7.26203472222222, 6.65711111111111, 6.40924305555556, 8.31249305555555,
  36. 8.32240277777778, 7.90968055555556, 8.42810416666667, 8.95072222222222,
  37. 9.70318055555556, 9.85254861111111, 11.4095486111111, 9.16079166666667,
  38. 11.0159375, 7.79172916666667, 9.49722222222222, 5.92420138888889,
  39. 7.548875, 5.80441666666667, 6.89536805555556, 5.63051388888889,
  40. 6.84440277777778, 6.24066666666667, 7.36372222222222, 5.56700694444444,
  41. 6.445375, 6.19682638888889, 6.48272222222222, 6.77991666666667,
  42. 7.06177777777778, 7.8626875, 8.15279861111111, 8.97194444444444,
  43. 9.16716666666667, 10.5453472222222, 10.578875, 10.2557361111111,
  44. 11.9828888888889, 10.7251736111111, 12.4082430555556, 11.4249166666667,
  45. 12.098875, 11.8695, 12.2577361111111, 11.9151319444444, 12.2548541666667,
  46. 11.4099513888889, 11.7989097222222, 11.0229027777778, 11.2145,
  47. 10.8331736111111, 11.2092083333333, 11.2498819444444, 11.6196180555556,
  48. 11.2190416666667, 12.3128958333333, 9.60270138888889, 11.6026111111111,
  49. 8.4145, 9.01735416666667, 9.02104861111111, 9.877, 9.37068055555556,
  50. 10.2840416666667, 9.59797222222222, 10.5512777777778, 9.90951388888889,
  51. 10.4284722222222, 10.5020347222222, 11.006375, 11.0382777777778,
  52. 11.6005972222222, 11.0032222222222, 11.5620277777778, 11.9508472222222,
  53. 12.2476041666667, 11.3607430555556, 12.1189722222222, 10.9222222222222,
  54. 12.3047291666667, 9.87395138888889, 10.6235972222222, 10.3749027777778,
  55. 10.7676666666667, 10.7136388888889, 11.9661875, 9.41833333333333,
  56. 11.2963611111111, 9.05438194444444, 11.3133680555556, 9.29961805555556,
  57. 10.7892638888889, 9.13478472222222, 9.78157638888889, 10.0053055555556,
  58. 10.2844236111111, 11.2174722222222, 11.3485555555556, 10.5509652777778,
  59. 11.3064583333333, 10.9914722222222, 11.1644097222222, 11.0532986111111,
  60. 11.3837847222222, 10.6440347222222, 11.1492361111111, 10.2033819444444,
  61. 10.906625, 8.8553125, 11.3004097222222, 7.00258333333333, 9.425375,
  62. 7.39536805555556, 8.19988888888889, 6.81268055555556, 7.90986805555556,
  63. 6.98128472222222, 7.428125, 7.55517361111111, 8.1825625, 7.2924375,
  64. 7.61071527777778, 8.64228472222222, 8.76127083333333, 9.1655,
  65. 9.96297916666667, 9.12442361111111, 9.43466666666667, 9.36382638888889,
  66. 9.8573125, 7.73863888888889, 9.06475, 6.86245138888889, 7.268,
  67. 7.28360416666667, 8.03669444444444, 7.46970833333333, 7.80116666666667,
  68. 7.90136111111111, 8.77240972222222, 8.38903472222222, 8.69964583333333,
  69. 9.0513125, 9.15219444444444, 10.0108819444444, 9.94815277777778,
  70. 9.73478472222222, 10.1262916666667, 9.17533333333333, 9.6718125,
  71. 9.18674305555555, 10.1854375, 9.16936111111111, 10.0093680555556,
  72. 8.82407638888889, 9.55520833333333, 9.40145833333333, 9.53936111111111,
  73. 8.80591666666667, 9.36661805555556, 8.04295138888889, 8.44384722222222,
  74. 7.49107638888889, 9.37970138888889, 6.54023611111111, 7.47807638888889,
  75. 6.58035416666667, 7.26764583333333, 6.2289375, 6.55975694444444,
  76. 6.55352083333333, 6.78065277777778, 7.3918125, 7.81163194444444,
  77. 8.44518055555556, 9.05359722222222, 8.69172222222222, 9.11545138888889,
  78. 9.03126388888889, 9.05704166666667, 7.69626388888889, 8.15723611111111
  79. ), se = c(0.0684446311524084, 0.0467000338163894, 0.0637176833357706,
  80. 0.0776305990169131, 0.0613172362964014, 0.0736477668094532, 0.0594146388107844,
  81. 0.0431317487238723, 0.0968421264052432, 0.0419846933555372, 0.0323949155766505,
  82. 0.0403999999512505, 0.0426753087710917, 0.04431240712146, 0.0685862399045098,
  83. 0.0938783494125271, 0.0908408653823437, 0.0266915764154258, 0.0953725248690779,
  84. 0.105824891014135, 0.0990518606457859, 0.0410508152091538, 0.06043250261061,
  85. 0.0786330141878098, 0.0652480589576939, 0.0696827474348336, 0.0535061067197427,
  86. 0.0741190834496391, 0.0701916014695556, 0.0735485250204054, 0.0458221447421354,
  87. 0.0404458878036966, 0.0359588714516765, 0.0265982467327904, 0.0394296624495542,
  88. 0.0427714222970342, 0.0221140283854524, 0.0384965936341351, 0.0782688867075244,
  89. 0.0247393796918376, 0.0382087990790521, 0.0475342767536421, 0.0479178691306617,
  90. 0.0416523624760465, 0.0559318984651426, 0.114747800062252, 0.0320540477244992,
  91. 0.0779599317541411, 0.0221973209356683, 0.0636192257453038, 0.0292632222064691,
  92. 0.0333513179466007, 0.0355241640648164, 0.0432614608297386, 0.0385724506235991,
  93. 0.037197666261358, 0.0402237976528405, 0.0491513763461124, 0.0612527892079573,
  94. 0.0568085063216254, 0.0687175273844523, 0.0294338061276067, 0.0636253008809029,
  95. 0.0912451311678636, 0.0384309704673781, 0.0221135650232409, 0.0499603906142495,
  96. 0.0732352504170765, 0.0496670437447336, 0.0392916567790627, 0.0518190951599677,
  97. 0.0451864397798832, 0.0449763902345077, 0.0376188816790774, 0.0353480974964637,
  98. 0.043738970054825, 0.0261383200661106, 0.0233255404890622, 0.0552663268605627,
  99. 0.0244428852057977, 0.0366775615825884, 0.0262850373316128, 0.0611982483162885,
  100. 0.0555639414319732, 0.0795167431194783, 0.0484336444188232, 0.0580722515120101,
  101. 0.0511506568571288, 0.0695252481623487, 0.047652522678636, 0.131585829693245,
  102. 0.0672199989613813, 0.0873260603254596, 0.106620462572398, 0.0429877006186673,
  103. 0.097603870391481, 0.044142339472256, 0.0616361009000405, 0.0444095513606984,
  104. 0.0671187567217355, 0.0882926151247376, 0.0591708833159303, 0.0722583700480519,
  105. 0.0521225933874304, 0.0867660870781919, 0.0927982580298578, 0.0286604318185285,
  106. 0.0165244497961037, 0.0215911564780479, 0.0271993439034762, 0.033820745033179,
  107. 0.0287039604443949, 0.046389154669938, 0.030690306262797, 0.144168379023189,
  108. 0.0656152472038106, 0.0713610892589061, 0.118974858624633, 0.100850818065691,
  109. 0.0946383067351734, 0.0526435227437007, 0.105288190570795, 0.0590228559096837,
  110. 0.0641402829806298, 0.0292958864631491, 0.0528746278867272, 0.0201978238929515,
  111. 0.0140800818414335, 0.0442948108734875, 0.0693191947399626, 0.0144639720670529,
  112. 0.04125110776308, 0.0095055017186666, 0.0303011849671232, 0.0125464998735671,
  113. 0.0489897150266719, 0.0751016435113841, 0.0789040475999884, 0.0440936049937659,
  114. 0.0211301628036506, 0.0755702005272336, 0.0540554625729186, 0.073399126077332,
  115. 0.0956125027382699, 0.0414118860902094, 0.04826324636017, 0.0501105502933141,
  116. 0.0534612307413563, 0.0340665419463383, 0.0207860495534631, 0.0147775084791049,
  117. 0.0124121556684425, 0.0176000529860602, 0.0147256372851489, 0.0289876122297381,
  118. 0.0244796725385413, 0.0127740158302701, 0.0158948404114211, 0.0242120175935621,
  119. 0.0228040488763762, 0.015269392245137, 0.0334098420839866, 0.0147282950908125,
  120. 0.0352636487530586, 0.0524948183516093, 0.0327090521078623, 0.0427079881785599,
  121. 0.0573994892952229, 0.062267865006105, 0.0297565541281921, 0.0142227669754567,
  122. 0.0675999462075173, 0.00781516310476744, 0.0415315319925174,
  123. 0.0295436118614684, 0.0251498717608074, 0.0465980555908037, 0.0382177558832005,
  124. 0.0449598834488143, 0.0293496702138142, 0.0151548810889912, 0.0434966522028971,
  125. 0.0364841418068765, 0.0294334202340891, 0.0233746698992999, 0.00814461650841544,
  126. 0.0172561458993778, 0.0183854419086275)), class = "data.frame", row.names = c(NA,
  127. -188L), .Names = c("jday", "site", "temp.avg", "se"))
  128.  
  129. dTMP$site<-as.factor(dTMP$site)
  130.  
  131. plot(temp.avg~jday, data=dTMP[dTMP$site=='A',], type='p', col='blue',ylim=c(5,13), pch=20 )
  132. lines(temp.avg~jday, data=dTMP[dTMP$site=='A',], type='l', col='blue',ylim=c(5,13), lty=1 )
  133. with(dTMP[dTMP$site=='A',], polygon(x=c(jday, rev(jday)) , y= c(temp.avg+1.96*se,rev(temp.avg-1.96*se)), border=NA, col=rgb(0,0,0,0.2) ))
  134. par(new=T)
  135. plot(temp.avg~jday, data=dTMP[dTMP$site=='B',], type='p', col='red',ylim=c(5,13), pch=20 )
  136. lines(temp.avg~jday, data=dTMP[dTMP$site=='B',], type='l', col='red',ylim=c(5,13), lty=1 )
  137. with(dTMP[dTMP$site=='B',], polygon(x=c(jday, rev(jday)) , y= c(temp.avg+1.96*se,rev(temp.avg-1.96*se)), border=NA, col=rgb(0,0,0,0.2) ))
  138. legend(x='topleft',legend=c('A','B'), pch=c(20,20),col=c('blue','red'))
  139.  
  140. library (nlme)
  141. cs1 <- corARMA(c(0.2, 0.3), form = ~ 1 | site, p = 1, q = 1)
  142. m2a<-gls(temp.avg~site, data=dTMP, correlation=cs1)
  143.  
  144. summary(m2a)
  145. anova(m2a)
  146.  
  147. Generalized least squares fit by REML
  148. Model: temp ~ site
  149. Data: dTMP
  150. AIC BIC logLik
  151. 454.5677 470.6964 -222.2839
  152.  
  153. Correlation Structure: ARMA(1,1)
  154. Formula: ~1 | site
  155. Parameter estimate(s):
  156. Phi1 Theta1
  157. 0.8661159 0.2603911
  158.  
  159. Coefficients:
  160. Value Std.Error t-value p-value
  161. (Intercept) 8.550266 0.7183725 11.902274 0.0000
  162. siteB 0.684303 1.0159321 0.673572 0.5014
  163.  
  164. Correlation:
  165. (Intr)
  166. siteB -0.707
  167.  
  168. Standardized residuals:
  169. Min Q1 Med Q3 Max
  170. -1.6126448 -0.5530600 0.2409433 0.9841246 1.7444175
  171.  
  172. Residual standard error: 1.949408
  173. Degrees of freedom: 188 total; 186 residual
  174.  
  175. arima(dTMP_BmA, order = c(2,0,1))
  176. Call:
  177. arima(x = dTMP_BmA, order = c(2, 0, 1))
  178.  
  179. Coefficients:
  180. ar1 ar2 ma1 intercept
  181. 1.2600 -0.5131 -0.7269 0.7458
  182. s.e. 0.1694 0.0973 0.1769 0.0553
  183.  
  184. sigma^2 estimated as 0.2353: log likelihood = -65.62, aic = 139.24
Add Comment
Please, Sign In to add comment