What’s up in Helena?

A quick look at Montana’s current legislative session

Like so many others, I’ve struggled to find my political footing in the current national flux.  One goal I’ve tried to chase is to find ways to be more politically informed at the local level.  To that end, I spent last weekend at the Women’s Policy Leadership Institute, hosted by Montana Women Vote in Helena.  I was there to learn some specifics about what’s going on in the current legislative session, and the Institute’s leaders definitely delivered on that front. I learned about the budget shortfall, bills promoting human rights, health care, and immigration policy in Montana, and found out about a bunch of local organizations working together to promote progressive ideas throughout the state. My summary (complete with bill names, links, and phone numbers wherever possible) is below.

      1. The Budget
      2. Bills to Help People
      3. Issues Indian Country
      4. Some truly shocking things about Montana criminal justice
      5. Local Projects

If you’ve forgotten who your legislators are, you can look them up using this neat map.

*A disclaimer: I know very little about policy, and this really was an imersive experience for me. I may have details wrong in some places, but I’m putting this out quickly since some of these issues are live right now.  I will continue editing as I learn more. If you see errors, please contact me!


The Budget

(quick overview)
Montana is facing a major budget shortfall in this session, largely due to an overly optimistic projection of revenues from oil and gas and corporate income tax, and to a lesser extent individual income taxes.  The state typically retains about $300 million in slush funds to deal with rainy day emergencies (and just this sort of missed projection), and fortunately those funds meant the legislature did not have to come back and change the budget midway through the last biennium. However, they’re gone now, and we don’t have additional rainy day reserves. Therefore, the legislature and the governor are looking for ways to constain spending in the state, while also restocking that reserve fund.

(proposed solutions)

Montana’s state budget is on the order of $2 billion per biennium. The budget consists of balancing appropriations with revenue. Governor Bullock proposed in his budget to deal with the shortfall by making $109 million in cuts to the state budget, reducing new proposal expenditures by $77 million, and raising an additional $150 million in revenue (through progressive tax increases on individual incomes, specifically adding a new top tax bracket and increasing the tax rate for individuals in that bracket by ~2%, as well as an ensemble of sin taxes, mainly targeting tobacco products, wine, and medical marijuana).

The legislature, unsurprisingly, is not interested in raising revenue (except possibly through sin taxes). They instead propose taking the entire $109 million Bullock proposed in cuts, and finding an additional $40 million to cut from the state budget. On top of this, and this is the big problem, the legislature proposes cutting $115 million from the state’s General Fund.

Now, here’s the issue. Money for most programs funded through the general fund receives about a 2:1 federal match, so that every dollar Montana puts in is matched by $2 from the feds. So, if we cut our contributions by $115 million to those programs, we will lose access to $345 million total.  This equates to roughly half a billion dollars in cuts.

(proposed cuts)

The legislature proposes to make a lot of the cuts from the Department of Public Health and Human Services (DPHHS). Approximately $100 million would come from health care, child and family services (though this area is likely to fare a bit better because Bullock has some pet early childhood development projects), services for the disabled, and largely, from senior and long-term care, who are currently scheduled to take on ~$80 million in cuts.

Higher Ed is also projected to take a big hit, to the tune of $23 million (here’s the talking point: counteracting that would require an estimated 23% tuition hike).

The Department of Corrections (who is apparently extremely strapped for cash) is projected to a $3 million cut.

(next steps)

There are a few steps forward here. First, it does look like the revenue projections were a bit pessimistic, and revenue numbers are still being updated. The legislature is currently working with an old set of numbers. They could move to update to the new projections, which would reduce the amount to cut some (I’m not sure how much). It is not entirely clear that they will do this, because some of them are just so fiscally responsible. This is a point on which constituent feedback may be helpful.

The legislature’s current set of cuts doesn’t yet meet the state’s monetary gap, so there is some possibility that they could cut further. Keep an eye on this over the next couple weeks as the legislative sub-committees complete their work and move forward.


Bills to Help People

The legislature is trying to do some good work on behalf of marginalized Montanans. This is a quick list of what’s been proposed, and how to support (or oppose) these ideas.

  1. State Earned Income Tax Credit (EITC)

The federal earned income tax credit is regularly tauted as a critical means for getting money into the hands of folks with very limited resources. Both the Montana House and the Senate are considering state EITC bills that would expand on the EITC for Montanans. The Senate’s bill (SB156) proposes a program at approximately 3% of the federal rate; the House bill (HB391) proposes a program at 10% the federal rate. The budget-knowledgables at this workshop suggested either of these bills (but especially the House bill) would be really helpful to low-income Montanans.

To support the state EITC, you can call into the legislative office, (406) 444-4800, and leave a message of support for the entire House tax committee (leaving a message for a whole committee is a thing you can do, apparently).

2. Montana Human Rights Act amendment

The legislature is yet again trying to revise the Montana Human Rights Act to cover discrimination on the basis of sexual orientation, gender identity, and gender expression (these are not currently included in our state Human Rights Act). The proposed amendment is HB417, sponsored by Kelly McCarthy. The Montana Human Rights Network, and Montana Women Vote are both collecting stories of discrimination. If you have one that you’re willing to put forward, please contact them.

3. Paid family leave

The federal Paid Family Leave Act (FMLA) ensures unpaid family leave for employees in companies with over 50 employees. That means that 2/3 of Montanans are NOT covered by the FMLA. The state is trying to fill this gap with a 1% tax on employees (matched by a 1% tax on employers) that would allow everyone some paid leave. My understanding is that the planned payout is regressive, so that people at the lower end of the earning spectrum receive nearly 100% of their typical pay while on leave, whereas folks who earn more are paid out a lower rate (bottoming out around %60 of typical pay, I believe). Paid family leave is covered in HB392, which is currently with the House Business and Labor Committee. Call the legislature at (406) 444-4800 and leave a message for the House Business and Labor Committee in support of paid family leave, or provide stories to Time For Montana.


Issues for Indian Country

The Native American communities in Montana are getting caught in the same budget downfall as the rest of the state (apparently tribes’ budgets are tied to the state’s budget, as opposed to being a separate entity, which I had not previously realized). Montana Budget and Policy Center‘s Heather Cahoon noted four major projects that are likely to get axed due to the state budget shortfall:

      1. Indian Country Economic Development Program (which aims to support Native American entrepreneurs through small grants)
      2. Native Language Emersion programs
      3. Tribal College assistance programs (these are largely funds to support non-native locals who attend tribal colleges and are not currently funded. As an aside, non-native students account for between 3 and 30% of students at Montana’s seven tribal colleges)
      4. Indian Country Suicide Prevention Program (btw, if you missed this nytimes piece on suicide and Standing Rock, I recommend reading it)

Access to good medical care remains a crucial issue for Native Americans in Montana, particularly for folks not living on reservations with Indian Health Service facilities (Lillian Alvernez, a law student at University of Montana reported that IHS is pretty inaccessible for approximately 40% of Native Americans in Montana).


Truly shocking things about Montana criminal justice

One theme at the meeting was criminalization of poverty. I think this isn’t a new idea for most media savvy liberal types, but there were a few special features to the situation in Montana that merit special mention.

First, the distribution of the people in Montana’s correctional system are skewed heavily toward detention and pre-release facilities, which provide programming (or, “programming”) to help rehabilitate (“rehabilitate”) people transitioning from the criminal justice system back to normal life. Most of these facilities (5 of the 7 pre-release facilities) in Montana are privately owned, and subject to the totally normal capitalist pull of maintaining a user base. Their incentive structure is not necessarily to set people up so that they leave the system forever.

Our criminal justice system costs roughly $200 million annually (compare that to the state budget above), in a system where 8 of 10 adult inmates are non-violent.

Prison poverty rates are absolutely unreal, with 37% of female inmates and 28% of male inmates earning less that $600/month prior to entering prison. 37%. Wow.

Yet incarceration rates in Montana are rising by about 200 sentences per year, even as sentence deferral rates are falling across the state.

There is currently a bill in the legislature (I missed the bill number:( ) to tie police academy funding to a court service fee. This is bad, since it financially incentivizes police to arrest people and have them processed by the courts.



Local Projects You Could Support

I learned a lot about different local efforts that I’d previously missed or overlooked. Here are a few.

Montana Women Vote

Montana Budget and Policy Center (they have great reports)

Montana ACLU (underfunded, understaffed, doing great work to support local action)

Gallatin Progressive Action Network (for your Bozeman-area protest-related needs; find them on facebook)



Montana Immigrant Justice Alliance



Montana Racial Equity Project (based in Bozeman; fighting racism in MT)

Montana Human Rights Network (doing good things for Human Rights in Montana)


Choice stuff

Planned Parenthood of Montana

Susan Wickland Fund (provides grants to help women in mountain west get to providers; based in Livingston)


Advancing Native American Communities:

Western Native Voice (based in Billings; supporting Native American communities)

Indian People’s Action

Native Generational Change

ACLU of North Dakota (has way too few people — what I heard was two people total — on the ground right now…)



On this “morning after”, let’s think about Roe.

Last night was surprising, shall we say, for this college-educated white woman. Call me a naive fool, but even as a lifelong resident of Mountain Time, I thought today we’d be looking at an historic “first” of a different kind. And yet, here we are. There’s nothing to do but move forward and try to get onto the same page as a society. In keeping with that theme, I’m writing about a topic I’ve considered writing on in the past, only to back away from it as too controversial. Yet this is the issue I thought about first while Florida, Ohio, North Carolina, and Michigan were turning from purple to pink. On this morning after, let’s think about Roe.

Roe, as we’re regularly told by the political left, is important for the health and safety of the mother and child. Less-discussed is another extremely important and perhaps further-reaching role: its critical place in the self-actualization of women. It is that somewhat bastardized side of Roe, which makes us all just a little squeamish, that I’ll talk about here.

Like many women, I spent the decade after college in training*, first as a statistician, and then as a disease ecologist. At 33, I’m currently one month into my first “real” job. I love my work, but that work hinges on my ability to occasionally work long, strange hours, and will likely one day require me to drop everything and move. Extreme geographic and temporal flexibility, after all, is essentially a professional expectation of early career practitioners, be they in academia, law, business, medicine, research, or any number of other professional domains.

Taking on a decade of education at serfian wages and under uncontrollable scheduling is a gamble with a relatively stable payout, assuming one can accrue enough work-years and appropriate advancement after finishing. Graduate and medical students (and, to a lesser extent, law students and MBAs) essentially make an up-front payment, on the bet that they will have sufficient time to reap the rewards. For female students, however, pregnancy risk alters the scope of that bet: reproduction, career advancement1,2, and financial stability3 are inexorably (and negatively) linked, due to the cost of delivering and rearing a child, plus –importantly – the lost financial opportunity associated with being a mom (expenses and educational commitments, in fact, are commonly cited as determining factors by women who actually have abortions4,5). As a consequence, many early-career women are effectively betting that we can make it through graduate school, plus put in enough post-grad time to justify the cost, prior to getting pregnant. And unless the Republican platform was a ship of lies (…), we’re unlikely to decouple pregnancy and women’s long-term financial stability any time soon.

Men, as far as I can tell, are largely blind to the pregnancy-finance risk conundrum facing their female friends and colleagues. I suspect it’s pretty easy for men (even well-intentioned, intelligent men) to remain blissfully naive6 to a very present worry for many young women: if a man isn’t aware of his partner’s cycle, he really has no means of counting minutes, waiting, and worrying anyway. Controlling reproduction (or rather, not-reproduction) remains very much women’s work, and it is work that women largely take on alone.

Women’s pregnancy-finance risk conundrum is compounded by evidence that women are more risk-averse than men to begin with, and are much more prone toward adopting strategies to hedge. It seems to me that risk hedging in the pregnancy-finance conundrum takes three obvious forms. One is to not have sex (which is a nice idea, but somewhat impractical given the increasing age at marriage for professional women, paired with increasing time until most couples WANT their first child); the second is to religiously (ha ha) use birth control (which most women do, and which is actually more of a pain then I think most men realize, but see 7); the third is to accumulate financial resources as quickly as possible so that the costs of an unexpected child could be borne (with the trade-off that jobs that pay well initially often do not allow for the skill-development necessary for upward mobility).

Let’s put aside Option One as unrealistic. Academia and modern society both implicitly and explicitly frown on Option Three, with (I suspect) the consequence that many, many women opt for Option Two. And Option Two, of course, isn’t fool-proof.

Enter Roe. Roe is the silver bullet, an essentially fail-proof method for dodging the pregnancy-finance spiral (until you’re ready for it). In its absence, I suspect more women would feel compelled to choose early nest-egg building over a long-term investment in education and career development. While it’s probably not true that no women want abortions, I suspect the vast majority of women regard them with the same fear and sadness that drives such emotional responses from the political right. That is probably the precise reason why abortions are increasingly rare. But having the option is imperative, since it gives women the freedom to invest in our futures along comparable axes to those experienced by men. As such, at this time, I believe Roe remains an essential tenet for the advancement and self-actualization of American women. I hope our leaders, new and old, can muster the empathy to see Roe in this light.

1. Mason MA & Goulden M (2004). Marriage and baby blues: redefining gender equity in the Academy. Annals of the American Academy of Political and Social Sciences 596: 86-103.

2. Wolfinger NH, Mason MA & Goulden M (2008). Problems in the pipeline: Gender, marriage, and fertility in the Ivory Tower. The Journal of Higher Education 79(4): 388-405.

3. Leung ML, Groes F & Santaeulalia-Llopis R (2016). The relationship between age at first bith and mother’s lifetime earnings: Evidence from Danish data. PLoS One 11(1): e0146989.

4. Finer LB, Frohwirth LF, Dauphinee LA, Singh S & Moore AM (2005). Reasons US women have abortions: Quantitative and qualitative perspectives. Perspectives on Sexual and Reproductive Health 37(3): 110-118.

5. Broen AN, Moum T, Bodtker AS & Ekeberg O (2005). Reasons for induced abortion and their relation to women’s emotional distress: a prospective, two-year follow-up study. General Hospital Psychiatry 27: 36-43.

6. Moore AM, Singh S & Bankole A (2011). Do women and men consider abortion as an alternative to contraception in the United States? An exploratory study. Global Public Health 6(S1): S25-S37.

7. Ricketts S, Klingler G & Schwalberg R (2014). Game change in Colorado: Widespread use of long-acting reversible contraceptives and rapid decline in births among young, low-income women. Perspectives on Sexual and Reproductive Health 46(3): 125-132.

* In the spirit of full disclosure, women’s studies and reproduction are not my areas of research. If I’ve misinterpreted the primary literature, or my arguments would be better served (or nullified) by other sources, please let me know.

Multi-scale models in infectious disease epidemiology

This post may be subject to revision, contingent upon on-going discussion. The current date of last modification is March 24, 2016.

Models that link within- and between-host processes are potentially important tools in disease ecology, but the disease research community remains somewhat divided about when and how they should be used. CIDD members discussed cross-scale modeling at a recent CIDD lunch. While the discussion was fairly free-wheeling and ad hoc, a few consistent themes emerged. This post draws on those themes, as well as a few of my own observations on the role – realized or potential – that multi-scale modeling plays in modern infectious disease epidemiology.

My background is in empirically-driven between-host modeling, and that baseline partially determines the challenges I identify here. This perspective is incomplete, and I hope this post will serve as a jump-off point for a broader discussion among researchers with more varied backgrounds.

The post consists of the following:

  1. Introduction
  2. What impedes multi-scale model development?
  3. When is a multi-scale model worth the effort?
  4. Parting thought

I’ve included a few references at the end, along with some comments from experts with other perspectives.

Introduction (next) (to top)

Multi-scale models use within-host dynamics to drive between-host epidemiological (or evolutionary) processes. This is implemented by treating state-transition rates in the between-host process (for example, transmission and disease-induced mortality rates in a classic SIR model) as functions of dynamics within the host. In their 2008 TREE paper, Mideo, Alizon, and Day note that multi-scale models are only essential in cases where reciprocal feedbacks link within- and between-host processes (for example, this occurs when SIR state-transitions depend on within-host dynamics, and also within-host dynamics depend on conditions at the population level). Mideo, Alizon, and Day observed that at the time they wrote, these reciprocal feedbacks were not actually included in the majority of publications on multi-scale models for pathogen evolution; in these cases, the resulting inferences could have emerged from examination of a single scale.

Yet in many situations, within- and between-host processes do have reciprocal feedbacks, whether our models account for them or not. For example, reciprocal dependencies exist in systems where within-host dynamics depend on initial dose, and where individual infectiousness or morbidity depend on both immune and pathogen dynamics. Since both of these dependencies are pervasive among infectious diseases, there is reason to think that multi-scale models may have broad applicability.

The apparent biological plausibility of multi-scale models raises two questions. 1) Which systems stand to benefit the most from multi-scale studies approaches? 2) If multi-scale models have so much utility, why aren’t we using them more frequently?

What impedes multi-scale model development? (next) (to top)

I see three clear barriers inhibiting multi-scale model development, each of which I describe below.

  1. Language and style
  2. Articulating model objectives
  3. Gathering and incorporating data

1. Language and style (next) (to top)

Within- and between-host modeling cultures are fairly distinct from one another, and for the inexperienced practitioner this poses a challenge for model design. Researchers trained in between-host dynamics come from a culture in which the vast majority of models derive from a single basic structure, Susceptible-Infected-Recovered (SIR). A rich knowledge of SIR-like models is crucial for disease ecologists, but comes with two important limitations.

First, SIR-focused modelers may not be particularly well-versed in other kinds of consumer-resource interactions, like mutualisms, competition, or predator-prey dynamics (for a good synthesis of consumer-resource models, see Laffery et al. Science 2015). This limits the set of model structures we draw from when considering the various co-stimulatory relationships governing interactions between pathogen and host immune responses.

Second, population-level disease modelers are used to working on models with relatively well-defined states. S, I, and R are good jump-off points for almost every between-host disease model. Although one state, pathogen population size, is an obvious fixture of most within-host models, constraining the host immune response to a limited number of appropriate compartments is really difficult. This issue is somewhat compounded by the immunology community’s emphasis on specificity: binning different groups of immune markers into single entities is counter to current research in that academic culture. I suspect that the juxtaposition of population-level disease modeling’s tendency to bin and immunology’s tendency to split states stymies many multi-scale modeling efforts before they even begin.

2. Articulating model objectives (next) (to top)

Disease researchers trained on between-host processes (or at least, the one writing this post) sometimes have a tendency to model first, and ask questions later. Because the equilibrium conditions of the SIR model and its various derivatives have been so extensively studied, just constructing and parameterizing an SIR-like model can offer a number of useful insights for a given system. Selecting which insight is most important from the outset of the modeling endeavor isn’t always imperative.

The ability to choose model goals after-the-fact is less available for within-host – and consequently, multi-scale – models. Modeling first and asking questions later isn’t feasible, since the set of possible model states and configurations in-host is so broad. Instead, within-host modelers tend to emphasize developing a specific question from the outset of an investigation, and selecting model states aimed at addressing that particular question.

3. Gathering and incorporating data (next) (to top)

There is a huge volume of published data on within-host dynamics for all kinds of pathogen-host combinations (I can attest that this is true from Mycoplasma ovipneumoniae, which very few people care about, so I’m willing to postulate that it’s probably true for whatever agent you’re studying right now). Making sense of these studies, however, isn’t trivial, and this is especially true for an outsider with limited microbiological and immunological literacy. Here’s why: unlike many between-host datasets, which are subject to relatively few (well-studied) assumptions with respect to data collection (e.g., sightability, population closure, etc.), the assumptions underlying many within-host datasets – which often directly manipulate the host’s health through gene knock-outs, particular nutritional regimens, and chronic stress – are daunting. Understanding and adjusting for potential biases also requires a reasonable working knowledge of a wide range of immunological and microbiological methods and their idiosyncrasies (how accurate IS that ELISA, anyway?), a good understanding of the study conditions, and an appreciation for how those conditions differ from conditions in nature. Furthermore, longitudinal datasets, or other data that capture temporal dynamics within hosts are relatively rare in immunology. For ecological modelers whose data gold standard is rich, replicated timeseries, a lack of analogous data within the host seems acutely problematic.

However, this need not be the case. Other data structures (especially multi-sectional datasets) do contain relevant (albeit sometimes less) information, and replication under slightly varied conditions is exactly the pretext for Bayesian inference. The empirical challenge – largely solved through partially observed Markov process modeling, approximate Bayesian computation and other approaches, still rarely implemented – is to figure out ways to appropriately leverage these less-than-ideal datasets. Conflicting attitudes about when particular datasets can and should be used may impede this process, but the statistical infrastructure for these models already exists.

When is a multi-scale model worth the effort? (next) (to top)

Despite the challenges, I am convinced that multi-scale models are crucial to advancing infectious disease epidemiology. Here’s why:

We’re currently very good at modeling epidemics for which Infected is Infected is Infected. Despite the clamour surrounding Ebola modeling, I think between-host models generally perform very well for acutely immunizing infectious, in which the preponderance of heterogeneity is attributable to host behavior.


With the notable exception of HIV (perhaps a special case, given that so much heterogeneity in secreted load can be explained by infection age), we do not have good population-level models of most chronic pathogens. Population-level epidemiology of pathogens like Herpes Simplex Virus, cryptosporidium, or tuberculosis that are typically latent but occasionally “flare up” (apparently at random) within particular hosts is not well-characterized*. Understanding transmission dynamics of these pathogens requires understanding the within-host processes that allow for re-emergence, and the population-level consequences of increased individual infectiousness. In short, what allows sporadic epidemic transmission of a pathogen that was locally endemic the whole time. The answer is sometimes pathogen evolution, but I posit that other factors may also contribute. Within-host complexities that lead to sustained infectious periods and pathogen re-invigoration has been well-reviewed (for example, see the excellent review by Virgin et al. 2009), but I have yet to see these ideas extensively translated into population-level predictions.

I suspect that multi-scale modeling may be most useful at a different point in a research program than between-host models. While between-host models often occupy the ends of epidemiological projects (e.g., though assessing R_0’s sensitivity to particular aspects of the disease process), and lead directly to management recommendations, multi-scale models might be best used to generate hypotheses, identify key pieces of missing information, and constrain the relevant space of unknown parameters for future experimental investigation. As a consequence, model outputs might shift from estimation-focused between-host models to a focus on mathematical sensitivity. This transition will require some flexibility on the part of practitioners; however, although the benefits are not entirely clear from the outset, multi-scale approaches hold a lot of promise, especially for pathogens that continue to defy clear population-level modeling.

Parting thought (next) (to top)

Perhaps the real challenge, then, is that multi-scale models might be most useful in situations that defy standard assumptions both within and between-hosts. I think that the following kinds of infections fall into this group:

  1. infections with relevant spatial compartmentalization within the host (e.g., HepC, HSV)
  2. epithelial infections (e.g., infections of the respiratory, GI, or reproductive tracts)
  3. infections that stimulate autoimmune elements of the immune system
  4. infections in which immune regulation misfires (either by over- or under-response)

References (next) (to top)

Mideo N, Alizon S, Day T. 2008. Linking within- and between-host dynamics in the evolutionary epidemiology of infectious diseases. Trends in Ecology and Evolution 23(9); 511-517.

Handel A, Rohani P. 2015. Crossing the scale from within-host infection dynamics to between-host transmission fitness: a discussion of current assumptions and knowledge. Philosophical Transactions of the Royal Society B 370; 20140302

Lafferty KD, DeLeo G, Briggs CJ, Dobson AP, Gross T, Karis AM. 2015. A general consumer-resource population model. Science 349(6250); 854-857.

Virgin HW, Wherry EH, Ahmed R. 2009. Redefining chronic viral infections. Cell 138(1); 30-50.

Comments and discussion (to top)

Jessica Conway made the following comment, which I’ll quote directly:

[You claim that] “Understanding transmission dynamics of these pathogens requires understanding the within-host processes that allow for re-emergence, and the population-level consequences of increased individual infectiousness.” It doesn’t really, does it? Take HSV. To understand population-level spread, you need outbreak pattern data, length, duration, and frequency. Do you need to know what drives them? Probably not to understand spread. Maybe yes, to determine interventions that minimize spread, but for that you pretty much want to stop flare-ups, and the between-host scale is not helpful. To understand evolution however likely requires both scales.

This is a really good point. In-host models are probably sufficient, so long as initial dose isn’t a critical determinant of infection outcome in the host. If dose is critical, then I think there is an opportunity for feedback across scales that could make multi-scale models helpful.

Field posts, 2015

Saturday, Sept 19th



Tuesday, May 12th

Transit day! I left Logan to keep the Hells Canyon project under control and headed back through the Lochsa to Bozeman this morning to prep for the Montana leg of this project. Before I left, though, Logan and I got to see 149.849 deliver her lamb on Bracken Point in Asotin.

Bracken Point is a pretty sweet place to view. The ewes like to tuck into some pockets in the cliffside low on Bracken Point (a place Miggy and I dubbed the “lamb caves” during the 2013 iteration of this project). We can watch the sheep at the lamb caves from a mesa about half a mile away. From that vantage point, we’re above them, so we get a really excellent view.

149.849 was at the lamb caves alone yesterday. In my mind, “alone-ness” and “lamb caves” bear a high association with lambing at this time of year, so we speculated she’d likely deliver soon. Logan went out at 6:30 this morning to check on her, and saw that her water had broken. We both went back together after that to watch. By 7:30, the ewe was having regular contractions on about 5-minute intervals. We could see lamb feet starting at 7:54, and the lamb was born exactly an hour later, at 8:54 this morning.

149.849 was up and down quite a lot while we watched. Her contractions were very clearly obvious, and she experienced most of them while bedded. Her tail and hindquarters would flex in undulating cycles of maybe 5 sets of 15 seconds.  After each round of contractions, she would stand and stretch.

Following birth, the lamb was moving and struggling almost immediately. It had its head raised within a minute of birth, and was standing (albeit not very well) when we left at 9:06. The ewe passed the placenta about 8 minutes after the birth, and ate it and umbilical cord pretty much immediately. She did not immediately pay a lot of attention to the lamb. While we almost always see ewes bedded against their lambs in the first few weeks, she did not immediately associate with it. I assume this was because she was still in pretty intense physical trauma from the birth, but I’m not entirely sure. Logan’s going to give me an update after a day or two.

Sunday, May 10th

We camped at Mountain View last night. No new lambs there, but we did see a very pregnant-looking 1160 off by herself. I’m hoping Logan will find her with a lamb in a couple days! We spent a while watching 779, 739, and an uncollared ewe with their lambs across from the Grande Ronde Lodge. The lambs are SO small. They barely look like sheep right now. 739’s in particular seemed to have quite a lot of energy: it was climbing all over its mom.

We’re having a little trouble with the set of ewes on the big feature across from the mouth of Grouse Creek. They like to tuck back into a canyon that we can’t really see from anywhere. So far, we’ve been lucky and caught them when they’re on the face, but I’m not sure how long our luck will hold.

Friday, May 8th

New lambs for 779 and 1409 in Mountain View!

Thursday, May 7th

New lamb for 145 in Asotin today (the first Asotin lamb this year). She was just above the Y along Asotin Creek Road. We got a look at them when the lamb was very new: dry, but not moving very well yet. Saw everybody else in Asotin, but no other lambs. Logan is headed to Mountain View alone tomorrow for  a little baptism by fire (and I get to spend the day cleaning data…).

Wednesday, May 6th

First lamb! 150.739 had a lamb at Mountain View today. I know she didn’t have it on Monday when I saw her, since she was the sheep out for the spontaneous swim in the Grande Ronde. I suspect she had it Monday night or yesterday, since it was up and moving (and dry) when we saw it this afternoon. It is SO SMALL: it can stand under her fully, with its head up, and not hit her belly.

We also saw everybody else at Mountain View, except 151.409 (who’s probably off lambing somewhere).

We got to see the house at the mouth of the Wenaha River that we can use if the weather is too wet to camp. It’s pretty sweet: gas stove, gas fridge, gas lanterns (that we won’t use, since they’ve burnt the paint on the ceiling…). It also has a really neat water-heating design: there’s a tank that afixes to the side of the wood stove with a spigot on it. You can fill the tank with water and then light the stove, and the stove will heat the water.

The Oregon refuge manager told me that the Grande Ronde is currently running around 2000 cus. Typically at this time of year, it’s 5-6K. Could be a smokey summer.

Tuesday, May 5th

Today was Logan’s first day. We worked Asotin, and saw pretty much everybody (but it did take us the whole day). There were 17 sheep on Bracken Point (including the ram 1240), six down with 145 at the forks of Asotin Creek, five in Charley Creek with 1389, and nine up at Lauffer Spring with 1851 and 1220. The Lauffer group took us forever to find: they were tucked back into some draws, and we could occasionally hear 1220, but weren’t getting 1851 or 360 (who was also with them).

Monday, May 4th

Today (my last solo day — my tech Logan starts up tomorrow), I worked Mountain View and Shumaker in Black Butte (see the map). Saw almost everybody in both places. I got lucky with this group of five yearlings and two ewes in Mountain View: they had no collars among them, but they were right in the middle of the road (looking for trouble).


I also saw this little guy cruising along on some cliffs just above them (it’s possible he was pushing the sheep group, although they didn’t seem particularly worried).


Still no lambs, anywhere, but who knows what tomorrow will bring!


Friday, May 1st

I worked Mountain View today, and it took me the whole day to get all the collars. No lambs yet. I’m still trying to figure Mountain View out. There’s a great paved road along the Grande Ronde, which gives us pretty easy access from below… but watching sheep from below is a bad idea. They like to lie right along edges, so that they’re stacked in. You might see a couple horns from below, when from above you could see that animal, plus five more bedded behind it. There’s an option to get up high at Mountain View, and look across the Grande Ronde from Grouse Flats (see the map here), but today a big group of ewes were tucked back into a draw that wasn’t visible from across or from the river. I caught a glimpse of them from the south side of the Grande Ronde, but I was probably looking over several miles. I could see collar colors, and make a reasonable guess at demographic group, but that was about it. Fingers crossed that they don’t decide to rear their lambs there!

Thursday, April 30th

Still no lambs at Asotin, but I did miss 1851 today, and she went early last year. I checked the spot where she lambed last year, up Lick Creek, but didn’t get any signal. Everybody else was cooperative, which was good since I needed to work fast today: Frances and I met for a few hours in the late afternoon to work on her paper describing the introduction of the new Mycoplasma ovipneumoniae strain and the ensuing epidemic at the mouth of the Grande Ronde in Black Butte last summer, so my field time was a little abbreviated.

We talked briefly about the count discrepancies we’ve been getting at Shumaker. Apparently Frances and Nick saw eight rams with the lone Joseph Creek ewe on their flight this week, and they think some of those rams are from Shumaker. Joseph Creek is the third part of the Black Butte population (along with Shumaker and the mouth of the Grande Ronde). The ewes don’t mix a whole lot between the three groups, but they share rams, and apparently also share pathogens. Frances and Nick think maybe the ram lamb who got lost his horn in the capture had been with the Shumaker rams, but when the rams split (for whatever reason) for Joseph, the ram lamb stayed behind and switched to the ewe group. We really need a closer look at that little guy (and a picture!) to determine if he’s actually ear tagged. All three of us thought we saw one, but none of us are totally convinced since the light was pretty bad (and I’ve only ever seen them from pretty far away).

Frances thinks Shumaker’s going to be rough to work during lambing. The ewes have cleared way out of Cassady’s land to lamb in the past, moving up river of Shneider Bar, or down further onto private land where access is pretty bad. She suggested we have Shumaker be our third priority, after Asotin and Mountain View (which she thinks should be tied for first priority; I’m not 100% convinced, since I like Asotin so much. I’d put it first, but only slightly, over Mountain View).

Frances also asked about 1891. I guess she didn’t move for a few weeks in March, and Nick had noticed that she was pretty thin (I’ve also noticed this, and her molt is ahead of everybody else’s). Frances and Nick thought she might have miscarried then. We took a look back at her old records (she’s been collared since 2008, but her frequency changed from 431 to 1891 during the winter of 13/14). She has a pretty bad track-record with lambs: none from 2008-2011, then one she lost in 2012, then two successful rearings while I’ve been here.

I’m trying a new schedule tomorrow, starting at Mountain View, then Shumaker, and then Asotin (namely, 1851). We’ll see how it goes.


Wednesday, April 29th

I was in town this morning for a check-in call with Pete, my advisor, and spent about three hours there moving and cleaning data, circulating plots, and checking in with other collaborators. I went to Shumaker and Mountain View in the afternoon, and did okay: saw everybody at Shumaker (no lambs yet… but I am still seeing 15 sheep, 9 ewes, 4 yearling rams, 2 yearling ewes), and got all but two at Mountain View (no lambs there either). Mountain View is going to be tough. It took me two hours to get home from the Grande Ronde Lodge to Asotin, and that alone will make for long days. Planning to start camping there every other day or every third day to cut down on the drive.

Saw this pretty guy on my way home (I took a close-up of his tail before getting out of the truck, just to check).


On my trip home, I saw a group of eight bull elk (big guys!). I don’t often see them in this phase. Their antlers remind me of well-trimmed Bonzai trees right now. The pictures are pretty granulated because it was pretty dark when I saw them, but you get the idea :)


Tuesday, April 28th

A short day in Asotin today. I saw 38: 13 yearlings and 25 ewes. 444, 149.849 and the ram 151.240 were up on Bracken Point with a bunch of yearlings. I saw 145 and fifteen others at the Y, a group of 7 including 1891 up Lick Creek, and 1851 and four of her friends at the mouth of Dry Gulch. Lots of good follows and scans, but no lambs in Asotin yet.

Monday, April 27th

Spent today with Logan, Frances, and Nick visiting Mountain View and the field house at the Wenaha, and checking in on the ewes at Shumaker. We got everybody except for two tagged yearling lambs at Mountain View, and saw everybody at Shumaker. The Shumaker light was pretty bad while we were there, so we didn’t get 100% confirmation on age structure, but we definitely saw 15 sheep. Frances and Nick were really excited that one lamb that they thought might have died shortly after capture appears to have made it – they hadn’t seen him since the capture, but it looked like he was in with this group today. I’d been really confused about his identity: he’s missing half a horn, so I was calling him the “half-uni ewe” earlier this week, but they think he’s a ram lamb. I’m shooting for a better look next time I’m out there.

Saturday April 25, 2015

Today was a classic spring day: rain, snow, and hail, plus plenty of sun.  I saw all the Asotin ewes (except 1180…) this morning. A little bad news from those girls, though: 150.319 and 151.361 were both coughing (319 quite a lot, just a couple singleton coughs from 151.361). I’m pretty worried about how the Asotin lambs will do this summer, given the symptoms I’m seeing.  Fingers crossed that they pull through okay.

I ran across this group of elk up on Smoothing Iron at about 7:00 this morning.  They weren’t too bothered by me.


Saw my first chicks this morning too, at the park in Lewiston, and heard through the grapevine that some of the other Hells Canyon sheep populations have lambs on the ground — I guess spring is really here!



Wednesday April 22, 2015

First day at Shumaker today, checking in on the Black Butte ewes. Saw 15. They seem to have had reasonable overwinter lamb survival, but some of the ewes are pretty scrawny.

Added bonus: the drive from the Asotin field house on Smoothing Iron to Shumaker is absolutely gorgeous. Should make for fabulous morning commutes:)

Saw this beautiful pair of herons on the Grande Ronde this morning while sheep-watching.



Tuesday April 21, 2015

Yesterday, I did some preliminary relocations, but today I really hit the ground running on scan samples (9) and focal follows (8) in Asotin. Tomorrow, I head to Shumaker (trying out the back-roads through Anatone; we’ll see how it goes).

I started the day with a group of 6 ewes, 151.561, 151.891, 150.145, 150.319, 151.600, and 151.361, up Charley Creek. They were headed down from the top of Charley, and I watched them from across the drainage on the west side. I got good locations, molts, and body conditions for all of them, and spent a little over an hour with them doing scans and follows.

Next was the group of 16 on Asotin Creek Road,  just downstream of the mouth of Charley Creek (initially, this group was only a couple hundred meters from the other group, but they pretty definitively split when I started collecting data on the 151.561 group). This group had a bunch of yearlings, (at least) two two-year-olds (orange 15 and yellow 41), and 151.251, 151.180, 151.220, 150.360, 151.851 and 150.020. I had a rough time getting molts and body conditions on these guys and reading eartag numbers since they were so high. Tried again in the evening, but they’d moved into Dry Gulch and the light was awful. Got some reasonable scans and follows during my morning visit, however.

Spent the afternoon with the group of 15 on Bracken Point, including 151.240, 150.444, 151.915, 151.389, 149.849, 151.811, and a BUNCH of 2-year-olds and yearlings. 151.240 (a two-year-old ram) was coughing, so I spent a long time watching for symptoms here. I contacted Frances, who reminded me that 151.240’s tested positive for Movi in the past. We’ll be watching for more symptoms down the road. Fingers crossed!  Got some scans and follows as well.

I saw leadership today from three sheep in particular in Asotin: 151.561 (leading a group with 145, 1891, 1600, 1361, 319), 151.251 (leading a big group with 1220, yellow41, orange15, 360, 151.851, 1180, others), and 149.849 (leading a group with ram 1240, 1389, 1915, 1811, 444 and a bunch of yearlings).

I’m having trouble distinguishing the yearlings from the two-year-olds, especially when trying to sort out the eartags. Hoping this becomes clearer with practice (and molt), and maybe some two-year-olds will have lambs….


Monday April 20, 2015

First day! Found all Asotin collars except 444 in three groups: downstream end of Bracken Point (15), just upstream of mouth of Charley Creek on Asotin Creek Road (20), and just downstream of mouth of Charley on Asotin Creek Road (3; different time of day than group upstream from mouth of Charley). Group of three was 020, 1851, and Orange 15. Hit the “sort” button on the relocations sheet and messed up all the points, but I think I’ve got them straightened out now.

Charley group had 5 years and a couple 2-year-olds. Got some pictures of them playing in the middle of Asotin Creek Road (just like the ones near the end of last summer when they were lambs).


Identifying animals when there are so many collars will definitely be a challenge. It’s possible to read the ear tags at close range, but takes time. 360’s and 1561’s collar both have almost no tape left. 180 was in with the Charley group, but her collar isn’t transmitting.

Had one class I ram (151.240) in with the group on Bracken Point.

One set of scans and three follows. Scanning with IDs will be rough since a lot of the taped collars are losing their markings; hoping for better tomorrow.

Pair of bluebirds are trying to nest in the upstairs of the Schlee house (flying repeatedly into upstairs windows), and a flicker is scoping out nesting in the corner of the porch. Saw some swallows (probably barn but I wasn’t 100% sure). Seeing lots of quail, chukkars, and hearing lots of meadowlarks (but haven’t actually seen one yet). Couple nice tom turkeys driving around yesterday. About 95 elk in yard of Schlee house when I pulled up Sunday night. Group of 7 beautiful mule deer up above the overlook on my run this morning.

Approximate Bayesian Computation for Ecologists

ABC (approximate Bayesian computation) is getting a lot of hype in the blogosphere and, increasingly, in the ecological literature (e.g., Hartig et al., 2011; Csillery et al., 2010). The idea is pretty appealing: link your super-complex simulation explicitly to your data with something resembling (or not, as the case may be) statistical rigor.

I like ABC. It has the flexiblity to give me access to complex process models without having to go through the pain (or impossibility) of writing out my likelihood. But I see the warts: an approximate posterior is not ACTUALLY a posterior; we are not ACTUALLY doing Bayesian inference in the classical sense. And there are a lot of ways in which we could mess things up.

The goal of this post is to introduce ABC, go through a few of the ABC caveats, and talk about how to assess whether (or to what extent) assumptions are being violated.

We get ABC courtesy of population genetics (thanks Pritchard!), where it’s been under steady development since the late 1990s. Tina Toni’s Journal of the Royal Society Interface paper in 2009 brought it (and specifically sequential Monte Carlo approaches for ABC) into the biological modeling mainstream. There are a few references to approximate Bayesian computation in the statistical literature beginning in about 2005, but that community realy became interested in ABC between 2008 and 2010, due largely to a few high-impact crossover papers that focused on genetic problems (see for example Cornuet et al., 2010). It’s been on the rise in cross-disciplinary modeling projects since then.

Intro: What’s ABC?
ABC often gets applied in cases where data are generated via some complicated process. The complication could lie in either the biological process, the observational process, or both; doesn’t matter. The data generating model is so complicated that we can’t (or won’t) write down an actual likelihood. It’s assumed that despite this complexity in the process, we have some nice data that we think are produced by the process. The goal is to use those data to estimate the parameters governing the process.

Case study:
I work with data that are likely generated by some deviant of an SIR (Susceptible – Infected – Recovered) model. SIR is a dynamical system where individuals transition from compartment to compartment. The model can sometimes have a whole bunch of parameters, and it’s fundamentally non-linear (the transition of individuals from S to I depends on how many individuals are in S and I at a given point in time). The data are (often) a batch of time-stamped disease occurrence observations. The challenge is to map the dynamical process to the empirical data and gain some meaningful insights into the process parameters.

The ABC approach to doing that follows these steps:
(1) Specify a summary measure on the data which you want to recapitulate in simulation.
(2) Simulate the process at a whole bunch of locations all over the possible parameter space.
(3) For each simulation, calculate a summary measure that is structurally identical to a summary measure calculated on the actual data. Reject the parameters that go with simulations whose summary measures don’t match the summary statistic from the actual data. The density of parameter values that weren’t rejected form an approximate posterior distribution for the parameters of interest.

There are a lot of flavors for implementing (2). The geneticists started off with a very simple rejection algorithm: simulate lots of parameters at randomly chosen locations, and throw out that ones that do a bad job, AKA ABC-rejection sampling. In the late 2000s, sequential Monte Carlo approaches came into vogue. These operate more or less like particle filters: (1) draw a batch of samples, (2) see what looks good, (3) draw more samples in the neighborhood of the samples that best matched the data in step 1, rinse, and repeat. Sampling at increasingly higher resolutions has the consequence of honing our examination of the posterior space to the regions with the highest approximate likelihood. Recently, methods are beginning to focus on incorporating hierarchical structures into ABC. By and large, these methods are still under development.

Issue number 1: Sufficient statistics are the best! (but what’s a sufficient statistic, again…?)

Several aspects of ABC (and especially that ever-desirable attribute, model selection) require sufficient statistics. In general, a sufficient statistic (for a given parameter) is a measure calculated on the data such that all relevant information in the dataset about that particular parameter is captured in the statistic. There’s an easy statistical trick to establish whether a given statistic is sufficient, called the factorization theorem, that allows a researcher to test sufficiency of any given statistic. However (and lucky for us), sufficient statistics for most commonly-used parameters are well-documented (for example, see the wiki list here). For a normal distribution with known variance, the sample mean is a sufficient statistic for the parameteric mean; for a Poisson distribution, the sum of the observations is a sufficient statistic for the rate parameter. For a sequence of independent and identically distributed (e.g., idd) Bernoulli trials, the sum of successes is a sufficient statistic for the success probability (given a known number of trials). The point is, we generally work with sufficient statitistics anyway, so running an ABC so that it matches sufficient statistics isn’t overwhelmingly difficult.

Issue number 2: How to explore your giant parameter space?
This is an area that has received a lot of attention. Toni et al. review this pretty thoroughly, but here’s the ten-second version. The basic strategy for many ABC approaches is to generate the posterior by running many very short chains, as opposed to a few very long chains. Sometimes (as in ABC rejection sampling), the chains are of length 1: draw a gazillion points from the joint prior, run a simulation for each draw, and build the posterior from simulations that produce statistics within some tolerance range of the empirical data. Sometimes (as in sequential Monte Carlo, AKA SMC), the chains are of length 3-5. Sample a bunch of initial points from the prior, simulate at those points, find the best matches to the empirical data, then resample from the prior in regions close to the best matches. Rinse and repeat a few times. Some new methods (e.g., Turner and Van Zandt, 2014) explore Gibbs-based approaches for sampling parameters in a hierarchical structure (where proposed parameters at different model levels respond to performance of parameters at other levels in the model). It’s worth noting that there are few (if any) examples of classical hierarchical modeling (as is implemented in lme4 in R, for example) meshed with ABC. I expect that will be an area of active work in the next few years.

Big questions in disease ecology??

Lately, I’ve been on a quest for big, open questions in disease ecology. I have lots of specific questions on my system, and some general developmental questions about the particular techniques I use, but fewer clear ideas about the big avenues that disease research will (or should) follow in the next ten years. In light of this apparent failing in creativity, I got the CIDD grad student group at Penn State to do some brainstorming with me a few weeks ago. Here are some of the questions we hit on, and in the last paragraph, some notable shortcomings.

Parasites as elements of ecosystems
1. What happens to vacated niches?
2. When does the dilution effect apply?

Questions about pathogen/host evolution
3. How do we measure “tolerance”, and can we rigorously describe selective drivers that control magnitude of host inflammatory responses?
4. What makes an evolution-proof vaccine?
5. What conditions lead to the virulence/tranmission trade-off, and what conditions produce evolution toward avirulence?
6. How do parasites transition to mutualists?

Questions about transmission within host populations
7. In general, what makes a super-spreader (and what counts more, super-spreaders or super-shedders)?
8. What selective pressures to pathogens impose on group size, and what’s our evidential basis for claims about that relationships?

Questions about within-host processes
9. To what extent is immune function an emergent property of the host’s microbiome?
10. Are there general features that describe how coinfection alters host immune response?
11. How does stress interact with immune response?
12. To what extent can we assume that within-host pathogen dynamics follow density-dependent growth processes? What are the costs of that assumption?

Questions about spillover
13. How does density of each host species shape spillover? Do general rules about density- and frequency-dependent transmission apply to spillover, and if not, what revisions are required?
14. How do we simultaneously manage for spillover and endemic disease in each host species (especially social hosts)?

Nina Wale pointed out that many of these questions are ecological in focus (she’s right: ecology is where I come from), and probably need to be reframed in an evolutionary light. I originally claimed that I was in seach of hypotheses, and Megan Greischar pointed out that modelers tend to spend lots of time aimed specifically at generating hypotheses. She has a point… but it seems to me like most models are built with some idea as to what hypotheses might eventually be relevant. In general, we found that asking big, general questions is uncomfortable. We much prefer asking specific questions about our systems. This is probably good: we should focus on the task at hand. On the other hand, I tend to get very focussed on the details of my work. As such, it’s important for me to occassionally try to nest my work in a broader scientific context. On the whole, it was an interesting discussion, albeit one that probably deserves quite a lot more of my time.

Planning Simulation Code in R

Planning simulation code

When I’m presented with a problem that I think I can solve by simulation, I tend to experience strong, adrenaline-laden urges to jump right into coding, and leave all the details to sort themselves out on the fly. Unfortunately, this is not a good idea (at least for me). The time saved by an adequate simulation plan definitely justifies at least identifying the simulation objectives, inputs, and outputs explicitly up-front. Since I write primarily in R, all but my harriest simulation projects happen there. I’ve received relatively little advice about how to structure simulations. This is my own protocol, so I’ll take full blame for all its anomales. It’s worked for me for a while.

Two minutes of planning
In a simulation intended to investigate the role age-structure plays in epidemic processes for a particular pathogen, I might want to look at maximum epidemic size under variable age structures, transmission rates, infectious periods, and disease-induced mortality rates. In this case, my overall inputs would be the conditions I want to vary, plus any other structural constraints, like an observed social contact structure, or population growth rate, etc., that I want fixed throughout all simulations. The output would be maximum epidemic size. I want to write my functions so that they can take my input values as arguments, and return my desired outputs.

An overarching strategy
My sim projects tend to (sometime inadvertently) follow an approximate Bayesian computational approach. In the simplest case, I simulate a dataset under some process — for me, it’s usually a disease invasion — with varying conditions (i.e., parameters). Next, I describe that particular simulated dataset with some predetermined summary statistics,  and finally I compare the simulated statistics to the same statistic calculated using field data. I throw out parameters that produce simulations with statistics that don’t match the observed data sufficiently well.

I break this procedure down into as many bite-sized chunks as possible. This usually gives me at least the following chunks (each embodied in its own function).

Function Utility
Priors.funDetermines parameter values that will govern a given simulationParameters for all prior distributionsRealized parameter values for a given simulation
DataSim.funSimulates one dataset (for me, usually longitudinal or timeseries) under a given set of parameters.Parameter values (user specified or drawn from prior)Either the entire simulated dataset, or a predetermined summary statistic
Batch.funSimulates replicate datasets under the same set of parameter values. (i.e., replicates DataSim.fun)Number of replications, parameter values to be passed to DataSim.funList or vector of length = number of replicates, containing summary measures for all replicate simulations
FullSim.funLinks Priors.fun, DataSim.fun, and Batch.fun into a single packageParameters describing the priors (to be passed to Priors.fun), Number of replications per set of priors (for Batch.fun)Matrix containing prior values and corresponding summary statistics

My preference is to write the DataSim function first.  This is typically the hardest bit. It contains the nuts and bolts of the biological process you’re exploring.  I find that I often identify overlooked aspects of the process as I write this function, so writing it first helps me clarify exactly how complex of a problem I’m tackling.  I often write smaller functions for critical and complicated steps in the data simulation (for me, these tend to be Disease Transmission, Survival, Recrudescence, etc.), so that I can deal with them in isolation. The Priors function can is probably the next-most challenging, since it involves all the philosophical and scientific growing pains associated with identifying reasonable priors. Batching over a set of parameters may not be necessary; I sometimes use it to quantify how variable my summary measures are under identifical parameters.  I often split FullSim.fun up into a few batches to pass to clusters, though this will depend on simulation size.

In general, I like to write all these functions in a file that contains the functions and nothing else (i.e., the “Source” file), and put everything else into another file (i.e., the “Run” file). More specs and suggestions for those in the next post. Below is pseudo-code for a source file (set up in R).

PriorsFun <- function(hyperpriors){
  priors.out <- fun(hyperpriors)
    # fun(hyperpriors) will likely use random number generators like runif, rnorm, etc.
    # the hyperpriors are the arguments that the random number generators need
  return(priors.out = priors.out)

DataSimFun <- function(priors.out){
  Data <-  
    #  might be complicated and include a loop through time
  summary.stat <- fun(Data)
  return(Data = Data, summary.stat = summary.stat)

BatchFun <- function(priors.out, n.reps){
  summary.stat.list <- vector("list", n.reps)
  # storage object for recording each sim's output
  for(i in 1:n.reps){
    summary.stat.list[[i]] <- DataSimFun(priors.out)$summary.stat
    # this could also be done with "replicate" but I like for-loops in examples,
    # since they're easier for me to understand
  return(summary.stat.list = summary.stat.list)

FullSimFun <- function(hyperpriors, n.reps, n.sims){
  param.list <- sim.out.list <- vector("list", n.sims)
  # storage objects for recording what parameters generated what summary measures
  for(i in 1:n.sims){
    param.list[[i]] <- PriorsFun(hyperpriors)
    sim.out.list[[i]] <- BatchFun(param.list[[i]], n.reps)
  return(sim.out.list = sim.out.list, param.list = param.list)


RStudio and R from PSU campus Windows computers

The motherboard on my laptop died a few days ago, and I’ve been using Penn State’s various campus computing labs to tide me through until my new machine arrives. Campus computers work fine searching articles, reading, etc., but figuring out how to comfortably run R (and RStudio) from them was somewhat more challenging. While the campus computers all have R installed, they don’t have RStudio, and I couldn’t even find the old-school default GUI that distributes with R.

Upon launching R locally on the campus Windows machines, a terminal (aka “Command Prompt” in Windows) opened with R running. This would have been fine in a pinch: since most of my code is functionalized, I could source it in and run things as I wanted in the terminal. However, I also wanted a decent text editor to comfortably edit the source files. I usually use vim on Linux (full disclosure: I really like RStudio, but running R from the terminal using the vim-R-plugin is my preference); on the PSU windows machines my text editor options were Notepad (no syntax highlighting), jEdit, and PSPad Editor (both found under All Programs -> Utilities -> Editors), neither of which had an R plugin for syntax highlighting available. Gah!

As luck would have it, I’ve been using the PSU HPC clusters for a few years, and attended a seminar on using R on the clusters earlier this fall. Usually when I use the clusters under Ubuntu, I ssh in from my machine, use scp to move a source file, a run file, and a PBS file that I’ve constructed and checked locally to my cluster space, and then run the PBS file on the cluster using qsub. I’d never tried editing code directly on the clusters themselves.

However, at the seminar I learned that one could access and interact with Hammer remotely: you can launch RStudio on Hammer, and edit and run code directly there, as I usually do from my own machine. AND this is easy to do from the public machines on campus! Thank goodness! Here are the steps:

1) Contact RCC and check to be sure you have space on Hammer (I think we’re all given space by default, but it’s worth checking just in case). In my experience, they’re really quick (minutes to hours) about setting up accounts when once requested.
2) Launch Exceed On Demand on the Windows Desktop. Here is a good description of using Exceed on Demand under Windows and OS X from RCC. You’ll be prompted to designate the following:

  • host (use hammer.rcc.psu.edu; to the best of my knowledge, Hammer is the only cluster for which this works)
  • User ID (your normal PSU ID)
  • Password (your normal PSU password)

3) You’ll be prompted to fill out a dialogue box about how you want to interface with Hammer (I’m currently using Xconfig = Seamless_Mode and Xstart = Terminal.xs); both seem to work fine.  Click “Run”.

4) A terminal will open on your machine.  To load and launch RStudio, do this:

module load rstudio

And lo, RStudio launches on your garden-variety PSU campus machine! No more procrastinating!