[R-sig-ME] Application of coxme

[davidcostantini at libero.it]

Dear All
I am trying to use the coxme function to compare survival among 6 experimental 
groups.
I have a fixed factor (group), a random factor (brood), and a censoring 
variable (0 = still alive;
1 = dead). My response variable is the age of the dead individual or of the 
individual still
alive at the end of the experiment.
I have tried to use this function:

coxme (Surv(time,status) ~ group, random=~1|nest, data=females)

Time would be the age, while status is the censoring variable.
R gives me an output and says that the random argument of coxme is 
depreciated.
I wonder if anyone of you may let me know if this function is correct and why
the random factor is depreciated. 

Cheers
David

>----Messaggio originale----
>Da: r-sig-mixed-models-request at r-project.org
>Data: 05/06/2013 12.00
>A: <r-sig-mixed-models at r-project.org>
>Ogg: R-sig-mixed-models Digest, Vol 78, Issue 11
>
>Send R-sig-mixed-models mailing list submissions to
>	r-sig-mixed-models at r-project.org
>
>To subscribe or unsubscribe via the World Wide Web, visit
>	https://stat.ethz.ch/mailman/listinfo/r-sig-mixed-models
>or, via email, send a message with subject or body 'help' to
>	r-sig-mixed-models-request at r-project.org
>
>You can reach the person managing the list at
>	r-sig-mixed-models-owner at r-project.org
>
>When replying, please edit your Subject line so it is more specific
>than "Re: Contents of R-sig-mixed-models digest..."
>
>
>Today's Topics:
>
>   1. Re: nlme: varIdent (Ben Bolker)
>   2. Re: Autocorrelation in a GAMM for nightly time series (Ben Bolker)
>   3. Re: Negative Variance (John Maindonald)
>
>
>----------------------------------------------------------------------
>
>Message: 1
>Date: Wed, 5 Jun 2013 01:14:37 +0000 (UTC)
>From: Ben Bolker <bbolker at gmail.com>
>To: r-sig-mixed-models at r-project.org
>Subject: Re: [R-sig-ME] nlme: varIdent
>Message-ID: <loom.20130605T031041-601 at post.gmane.org>
>Content-Type: text/plain; charset=us-ascii
>
>D Gill <oceancurrents at ...> writes:
>
>> 
>> Hello,  I am a student working on fisheries data using a mixed linear
>> model. I have a question about the varIdent function. I have information on
>> 140 fishers at nine sites, where the only random effect is for the site. As
>> I do not have equal variances between sites, is it correct to use the
>> varIdent argument for the sites? Do I lose any information by doing it this
>> way?
>> 
>> M1<-lme(resp~......,random = ~1|Site,weights = varIdent(form =~1|Site) ))
>> 
>
>  This seems reasonable to me.  However, note that you might
>want to consider (variable*site) interactions of any predictor variables
>that vary within sites (see e.g. Schielzeth and Forstmeier 2009) [e.g.
>if your measurement is total catch, effort is a covariate, and CPUE
>varies among sites, you would want effort|Site in your random
>effects specification].
>   I'm assuming you have a single measurement for each fisher -- 
>otherwise you should also add a random effect for fishers ...
>
>  Ben Bolker
>
>
>
>------------------------------
>
>Message: 2
>Date: Wed, 5 Jun 2013 01:35:51 +0000 (UTC)
>From: Ben Bolker <bbolker at gmail.com>
>To: r-sig-mixed-models at r-project.org
>Subject: Re: [R-sig-ME] Autocorrelation in a GAMM for nightly time
>	series
>Message-ID: <loom.20130605T032637-397 at post.gmane.org>
>Content-Type: text/plain; charset=us-ascii
>
>Andrew Digby <andrewdigby at ...> writes:
>
>> I'm analysing call counts of a nocturnal species in 
>> response to temporal and environmental effects. Over a
>> period of 3 years I divide each night into 10 equal-length blocks, 
>> and count the calls of each sex in each
>> block. 
>
> [snip]
>
>> gam(ncalls ~ sex + year + s(month) + s(time of night) +
>>  s(moon) + s(temperature) + s(rain) +
>> offset(blocklength), family=negbin(c(1,5)), data=dc)
> 
>> Where 'time of night' = 1-10 = the block number during the night,
>> and blocklength = the length of that block (this changes throughout
>> the year; hence the offset).
> 
>> Not surprisingly, the acf shows significant correlation between
>> adjacent time blocks, with a periodicity of 10. This is because
>> counts in a particular block will be similar to those in the blocks
>> immediately before and after, and also to the same block on other
>> nights (the species has a quite regular pattern of decreasing call
>> rates as the night progresses).
> 
>> To address these autocorrelation problems, I've tried various 
>> ARMA correlations, with correlation
>> between time of night blocks, grouped by day (or week) and sex:
>> 
>> gamm( .., correlation=corARMA(form=~TimeofNight | DayofYearSex), p, q) 
>>  (p=1:3, q=0:3)
>> gamm( .., correlation=corARMA(form=~TimeofNight | WeekofYearSex), p, q)
> 
>> This improves the ACF, but there are still correlations between
>> residuals at same time of night - e.g. peaks every 10 lags.
> 
>> 1) Can I rely on the ACF when my time series isn't actually
>> consecutive time bins, since it only includes each night (not full
>> days). This means, for example, that the time lag between
>> observations 9 & 10 (~1 hour) is much longer than that between 10 &
>> 11 (~12 hours = daytime)?
>
>> 2) (if yes to the above) How can I construct a correlation structure
>> in a GAMM which allows for correlation between adjacent time bins
>> within each night, and between the same time bin on different
>> nights?
>
>  I think you can use corCAR1, which is intended for continuous-time
>models, but it's limited -- it does only exponential correlations
>in time (analogous to corAR1). 
>
>  If you have enough data/computational power, you might try 
>a 2D spline (tensor product? I don't quite remember) to allow
>for a gradually changing nocturnal time pattern over the course
>of the year (i.e something like s(timeofnight,dayofyear)) --
>not entirely clear to me why you're fitting s(month) rather than
>s(dayofyear), unless you only have one observation per month?
>
>In general my experience is that when I find myself fitting
>high-order ARMA models (i.e. more than a couple of lags in total),
>it's because there's some systematic pattern that I failed to
>capture in the fixed-effects part of the model. Your mileage may
>of course vary.
>
>  Ben Bolker
>
>
>
>------------------------------
>
>Message: 3
>Date: Wed, 5 Jun 2013 11:45:13 +1000
>From: John Maindonald <john.maindonald at anu.edu.au>
>To: Ben Bolker <bbolker at gmail.com>
>Cc: r-sig-mixed-models at r-project.org
>Subject: Re: [R-sig-ME] Negative Variance
>Message-ID: <6C1A3D6B-E7C2-4136-A511-BFEE53D983CF at anu.edu.au>
>Content-Type: text/plain; charset=us-ascii
>
>A variance components model that has a variance structure 
>
>  block variance + plot (within block) variance + subplot (within plot) 
variance
>
>makes sense only if blocks take out some part of the variation, i.e., 
variation
>between plots within blocks is (in the absence of treatment effects) smaller
>than variation between plots in different blocks.  Similarly for subplots
>within/between plots.  
>
>If on the contrary, there is more variation between between plots within 
blocks
>than between plots in different blocks (this is likely to happen if there is 
a 
>nutrient or fertility or moisture gradient within blocks), then a model that 
has
>the form on the second line above will if allowed account for this by 
returning
>a negative block component of variance estimate.  It does this in order to 
get
>a plausible variance-covariance structure.
>
>Of course, once a gradient has been identified, it can be accommodated in the
>model.  This does not however undo all the malign effects of an unfortunate
>experimental design.
>
>John Maindonald             email: john.maindonald at anu.edu.au
>phone : +61 2 (6125)3473    fax  : +61 2(6125)5549
>Centre for Mathematics & Its Applications, Room 1194,
>John Dedman Mathematical Sciences Building (Building 27)
>Australian National University, Canberra ACT 0200.
>http://www.maths.anu.edu.au/~johnm
>
>On 05/06/2013, at 11:05 AM, Ben Bolker <bbolker at gmail.com> wrote:
>
>> John Maindonald <john.maindonald at ...> writes:
>> 
>>> 
>>> Negative variance estimates can be very useful in alerting that the
>>> variance-covariance structure is not what one expects.  Or they may
>>> allow the simplest way of specifying the overall variance-covariance
>>> structure, short of specifying the variance-covariance structure in
>>> some detail.
>>> 
>>> I was told of an experiment where the experimenters had chosen 
>>> blocks to be at right angles to the river bank, accordingly maximising
>>> between plot variance.  This came to light, in data analysed away
>>> from the scene of the original trial, when the block variance was
>>> estimated as negative -- a very useful diagnostic.  Certainly, one can 
>>> check on such a possibility by specifying a suitable variance-covariance 
>>> structure, but how many analysts will take that trouble?
>> 
>>  I don't quite get the geometry you're talking about, but
>> I take the general point that diagnostics are good and that
>> one wouldn't necessarily think to consider negative correlation.
>> 
>>> Or one has results from each of two eyes per person.  After allowing
>>> for any systematic left/right difference, are two eyes from the same
>>> individual more or less different than two eyes from different 
>>> individuals?  I doubt that there is a general answer that applies to all
>>> types of eye measurements.
>> 
>>  I find this one a little bit less convincing -- here it would
>> seem to be perfectly natural to fit a model that allowed for
>> positive or negative correlation.
>> 
>>  The fact remains that, whether or not it's a good idea,
>> this is very hard to do in nlme/lme4 for
>> structural reasons. Luckily people are suggesting alternative
>> packages.  (Anyone who would like to edit 
>> http://glmm.wikidot.com/pkg-comparison accordingly is welcome
>> to do so ...)
>> 
>>  I don't think "how do I estimate negative variances?" has quite
>> risen to the level of a FAQ yet, so I won't bother adding it
>> to http://glmm.wikidot.com/faq (although again, if anyone wants
>> to take the initiative to do so I wouldn't complain).
>> 
>>    Ben Bolker
>> 
>> _______________________________________________
>> R-sig-mixed-models at r-project.org mailing list
>> https://stat.ethz.ch/mailman/listinfo/r-sig-mixed-models
>
>
>
>------------------------------
>
>_______________________________________________
>R-sig-mixed-models mailing list
>R-sig-mixed-models at r-project.org
>https://stat.ethz.ch/mailman/listinfo/r-sig-mixed-models
>
>
>End of R-sig-mixed-models Digest, Vol 78, Issue 11
>**************************************************
>