I am relatively new to statistics and not formally trained but have been given a complex problem to solve and need some guidance. I realise that I am out of my depth a bit here but would appreciate whatever help I can get bearing in mind that there is no budget for this and as a result it is not possible to purchase software or hire consultants.
The Problem
The business I work for has a large number of mobile representatives that can be dispatched to a variety of different jobs. There are ~100 different job types and each job can be broken up into 4 different final outcomes. Each of these 400 outcomes requires an allocation of man hours to complete. I have a count of how many times each one of these outcomes occurred in each hourband for the past 5 years.
I have been asked to forecast how many of each outcome will occur in each hourband for the 28days from the present. The resulting forecast will be used to anticipate staffing requirements on an hour-by-hour basis. As a result the forecasts for each hourband need to been fairly accurate.
Factors
In my data there are clearly some yearly, weekly, and daily seasonal effects. In general each outcome is more likely to occur at certain times of the day on certain days of the week and with some yearly trends.
Each different outcome is likely to be related to the frequency of a number of different outcomes. i.e. if x happens then y and/or z are likely but a and/or b are not.
There are a large number of environmental factors that contribute to the frequency of each outcome. These can include, but are not limited to weather, sociopolitical, financial trends, one off events.
What I have tried
So far I have tried using simple auto.arima, holtwinters and ets forecasts. holtwinters ended up producing a flat line (i.e. 5 and hour for the next 672 hours). ets doesnt work because the seasons are longer than 24 intervals. auto.arima produced the best results but they were still a long way off being accurate.
It was then suggested that I try tbats() and provide it with multiple seasonal lengths. I achieved best results by giving it seasonal lengths of 8760 (1yr) and 168 (1wk). Frustratingly, these results are within 1% when viewed as a sum of all hourbands in a 1 month block but are anything up to 300% (avg 20%) off when considering each individual hourband.
Both of these approaches have been applied over an individual outcome rather than considering all possible outcomes (and their correlation to each other).
My thoughts so far
At this stage I feel like my two options are to either to find a way to use something similar to tbats() that will look at the relationships between the multiple different outcomes as well as the seasonality and forecast based on that information.
or
Abandon that approach for a Neural Network model. My understanding (limited) is that using the Neural Network approach I may be able to 'factor' for the multitude of unknown environmental factors without having to actually identify them. I know this is lazy but my feeling for the data is that there are going to be a fair few unknown factors to identify and forecasting them may end up being a job in itself (i.e. weather conditions)
The Question (finally)
What I am looking for is some guidance.
Considering the information above and the fact that I am pretty much limited to R, what is the best approach??
and
What are the basic steps I need to follow?
While I cant post my data online (due to my employers restrictions) I can send it an individual or two if someone was interested in giving us a hand to find a solution.