On May 20, 2011, at 3:54 PM, eeecon wrote:
> Hi,
>
> My code indicates there may be a bug in multilevel.
> I doubt this is actually the case, can anyone tell me what is wrong
> with my
> code?
>
> The data file for this code can be downloaded here:
> http://cameron.econ.ucdavis.edu/mmabook/mma15p4gev.asc
>
> Here is the code that generates the bug:
>
> ### rm(list = ls())
## Thanks a lot!
## That line should have been commented.
## To do otherwise says that the poster believes his problem is
## more important than anything else in the workspace.
>
> library("multilevel")
>
> #v1
> #fishing.data <- scan("mn.dat")
# This might be more readily reproducible with this line:
fishing.data <-
scan(url("http://cameron.econ.ucdavis.edu/mmabook/mma15p4gev.asc
") )
> fishing.data.array <- array(fishing.data, c(11, 4, 1182))
> fishing.data.array <- aperm(fishing.data.array, c(3,2,1))
> means <- aggregate(fishing.data.array[,,3:5],
> as.data.frame(cbind(fishing.data.array[,1:4,2])), mean)
>
> ## Table 1
> ## > means
> ## V1 V2 V3 V4 V1 V2 V3 V4
> V5 V6
> ## 1 1 0 0 0 35.69949 35.69949 97.80913 125.00316 0.2791948
> 0.2190015
> ## 2 0 1 0 0 30.57133 30.57133 82.42908 109.76335 0.2614444
> 0.2025348
> ## 3 0 0 1 0 137.52710 137.52710 41.60681 70.58409 0.2082868
> 0.1297646
> ## 4 0 0 0 1 120.64835 120.64835 44.56376 75.09694 0.2519077
> 0.1595341
> ## V7 V8 V9 V10 V11 V12
> ## 1 0.1593985 0.5176090 4.051617 4.051617 4.051617 4.051617
> ## 2 0.1501489 0.4980798 3.387172 3.387172 3.387172 3.387172
> ## 3 0.1775411 0.6539167 4.654107 4.654107 4.654107 4.654107
> ## 4 0.1771628 0.6914998 3.880900 3.880900 3.880900 3.880900
>
> #
> ## Note V1 and V2 above are the same. This is suspicious.
> ## In fact (1,V1) from Table 1 should agree with (1,V1) from Table 2
> below:
> #
It appears to me that you are having problems with your use of aperm.
>
> #v2
> fishing.data <- scan("mn.dat")
> dim(fishing.data) <- c(11,1182*4)
> fishing.data <- t(fishing.data)
> means <- aggregate(fishing.data[,3:5],
> as.data.frame(fishing.data[,2]*fishing.data[,6:8]), mean)
>
> ## Table 2
> ## > means
> ## V1 V2 V3 V1 V2 V3
> ## 1 0 0 0 95.85924 0.2717639 4.097600
Was not what I got.
> means
V1 V2 V3 V1
V2 V3
1 0.000000e+00 0.000000e+00 0.000000e+00 7.563181e-01 7.444081e-01
7.546115e-01
2 1.736111e-01 0.000000e+00 0.000000e+00 0.000000e+00 0.000000e+00
0.000000e+00
3 5.208334e-01 0.000000e+00 0.000000e+00 0.000000e+00 0.000000e+00
0.000000e+00
4 8.680555e-01 0.000000e+00 0.000000e+00 0.000000e+00 0.000000e+00
0.000000e+00
5 1.000000e+00 0.000000e+00 0.000000e+00 0.000000e+00 0.000000e+00
0.000000e+00
6 1.215278e+00 0.000000e+00 0.000000e+00 0.000000e+00 0.000000e+00
0.000000e+00
7 1.562500e+00 0.000000e+00 0.000000e+00 0.000000e+00 0.000000e+00
0.000000e+00
snipped about 1490 lines.
Would have been much more courteous to use:
head(means)
I think you should take the time to build a minimal example. Now that
my workspace has been trashed, I have lost enthusiasm for this project.
> ## 2 1 0 0 30.57133 0.2025348 3.387172
> ## 3 0 1 0 41.60681 0.1775411 4.654107
> ## 4 0 0 1 75.09694 0.6914998 3.880900
>
> --
> View this message in context:
http://r.789695.n4.nabble.com/multilevel-tp3539421p3539421.html
> Sent from the R help mailing list archive at Nabble.com.
>
> ______________________________________________
> R-help at r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide
http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.
David Winsemius, MD
West Hartford, CT