OxCal > Analysis > Examples

Analysis Examples


This part of the manual contains illustrative examples of input files.

Further examples are given throughout the Analysis section of the manual.


Multiple plots

Frequently, all that is required for presentation of a set of calibrated radiocarbon measurements is a plot of the probability distribution functions and ranges.

Such a plot is generated by performing multiple calibrations in one project file as outlined in the in getting started with multiplots. There are also various commands which can be used to help in the organisation of the plot. These are Axis(), Label() and Page(). The following shows a typical example.

 Plot()
 {
  Label("Area 1");
  Axis(-300, 300);
  R_Date("A",2023,20);
  R_Date("B",1961,20);
  R_Date("C",1999,20);
  R_Date("D",1966,20);
  Line( );
  Label("Area 2");
  R_Date("E",1954,20);
  R_Date("F",1936,20);
  R_Date("G",1948,20);
  R_Date("H",1925,20);
  Page( );
  Label("Area 3");
  Axis(-300, 300);
  R_Date("I",2044,20);
  R_Date("J",1912,20);
  R_Date("K",1730,20);
  R_Date("L",1820,20);
 };

See also [Getting started with multi-plots]


Combining radiocarbon dates

Radiocarbon dates that all relate to the same event (in other words where all of the samples should should all derive from the same radiocarbon reservoir at the same time) can be combined together before calibration. This is done using the R_Combine() function. The easiest way to enter dates for combination is to use the [Tools > Models > 14C date combination] tool. A couple of examples are shown here:

 R_Combine("Combination")
 {
  R_Date("A1",2023,20);
  R_Date("B1",1961,20);
  R_Date("C1",1999,20);
  R_Date("D1",1966,20);
 };
 R_Combine("Single year",8)
 {
  R_Date("A2",2023,20);
  R_Date("B2",1961,20);
  R_Date("C2",1999,20);
  R_Date("D2",1966,20);
 };

In this example you can also see how you can incorporate an additional element of uncertainty in the radiocarbon concentration. This might be useful for example in the case of single growth season material which will be expected to show some extra uncertainty (about 8 14C years) relative to the calibration curve which is usually based on decadal measurements (Stuiver et al. 1998).

See also [Tools | Radiocarbon calibration]


Plotting calibration curves

You can plot the default calibration curve by running an empty model. To plot specific curves you can use the Curve() command for the relevant curves. The Axis() command can also be used to specify the range of the plot. For example:

 Options()
 {
  Resolution=0.2;
  Curve="bomb21nh1.14c";
  RawData=TRUE;
 };
 Plot()
 {
  Curve("Bomb21NH2","bomb21nh2.14c");
  Curve("Bomb21NH3","bomb21nh3.14c");
  Curve("Bomb21SH3","bomb21sh3.14c");
  Curve("Bomb21SH12","bomb21sh12.14c");
  Curve("Kueppers04","kueppers04.14c");
  Axis(1950, 2000);
 };

will plot all of the bomb curves for comparison. The options set at the top of this will ensure that the binning resolution is 0.2 years and that the curve data-points are sent to the plotting utility. In the plot you can choose to plot these points (which shows where there in interpolation and binning) and you can also choose to plot against F14C rather than BP date.

See also [Radiocarbon calibration | Calibration Curves]


Finding the span of a group of events

In order to do this you need to specify that the dates are all part of a single phase. The easiest way to do this is to use the [Tools > Models > Phases] tool to create the phase. You can then enter the radiocarbon dates and use the Span() function to find the span of the dated events in that phase as in this example:

  Sequence()
  {
   Boundary("Start 1");
   Phase("1")
   {
    R_Date("", 3100, 20);
    R_Date("", 3000, 20);
    R_Date("", 3040, 20);
    R_Date("", 3080, 20);
    R_Date("", 3010, 20);
    R_Date("", 3030, 20);
    Span("Span of dates");    
    Interval("Duration");
   };
   Boundary("End 1");
  };

You can also use the Interval() query in the way shown here to calculate the interval between the start or the phase and its end (ie the duration of the phase). This duration includes events which have not been directly dated.

See also [Tools | Groupings | Queries]


Dealing with multiple phases

In many archaeological applications we need to deal with multiple phases. Such models are usually best entered using the [Tools > Models > Phases] tool. In this example we look at two different treatments of the same data. The measurements here are dates from hafting materials found attached to bronze axe heads (data from Needham et al. 1998).

In the first model (below left) we assume that the different phases are completely independent (overlapping phases) and simply try to estimate their start and end dates from the dates that we have for each type of material. In the second model (below, right) we assume that the phases are in the order deduced from the archaeological record (contiguous phases). In this second model we estimate the dates of transition from one phase to another.

 Phase("Bronze Age")
 {
  Sequence()
  {
   Axis( -3000, -400);
   Boundary( "Start of Acton & Taunton");
   Phase( )
   {
    R_Date( "OxA-5948", 3225, 65);
    R_Date( "OxA-4651", 3220, 80);
    R_Date( "OxA-5949", 3110, 50);
    R_Date( "OxA-6177", 3055, 50);
    R_Date( "OxA-5196", 3035, 40);
    Interval( "Duration of Acton & Taunton");
   };
   Boundary( "End of Acton & Taunton");
  };
  Sequence()
  {
   Boundary( "Start of Penard");
   Phase( )
   {
    R_Date( "OxA-5187", 3045, 55);
    R_Date( "OxA-1526", 3030, 100);
    R_Date( "OxA-5953", 3015, 45);
    R_Date( "OxA-5951", 2980, 45);
    R_Date( "OxA-4504", 2965, 65);
    R_Date( "OxA-5952", 2965, 45);
    R_Date( "OxA-5959", 2965, 45);
    R_Date( "OxA-5183", 2930, 40);
    R_Date( "OxA-4653", 2910, 55);
    R_Date( "OxA-5950", 2910, 45);
    R_Date( "OxA-5954", 3025, 55);
    R_Date( "HAR-2940", 3020, 70);
    Interval( "Duration of Penard");
   };
   Boundary( "End of Penard");
   Page( );
  };
  Sequence()
  {
   Boundary( "Start of Wilburton");
   Axis( -3000, -400);
   Phase( )
   {
    R_Date( "OxA-4656", 3005, 75);
    R_Date( "OxA-4503", 2930, 55);
    R_Date( "OxA-4502", 2925, 50);
    R_Date( "OxA-5036", 2920, 50);
    R_Date( "OxA-5197", 2910, 50);
    R_Date( "OxA-5955", 2900, 45);
    R_Date( "OxA-5035", 2900, 45);
    R_Date( "OxA-5034", 2890, 45);
    R_Date( "OxA-5956", 2850, 50);
    R_Date( "OxA-5198", 2820, 70);
    Interval( "Duration of Wilburton");
   };
   Boundary( "End of Wilburton");
   Page( );
  };
  Sequence()
  {
   Boundary( "Start of Ewart Park");
   Axis( -3000, -400);
   Phase( )
   {
    R_Date( "OxA-5957", 2810, 45);
    R_Date( "OxA-4716", 2780, 50);
    R_Date( "OxA-4654", 2765, 45);
    R_Date( "OxA-5976", 2740, 45);
    R_Date( "OxA-4652", 2720, 45);
    R_Date( "BM-798", 2704, 50);
    R_Date( "OxA-5962", 2685, 60);
    R_Date( "OxA-6176", 2655, 50);
    R_Date( "OxA-5977", 2620, 45);
    Interval( "Duration of Ewart Park");
   };
   Boundary( "End of Ewart Park");
   Axis( -3000, -400);
  };
 };
  Sequence("Bronze Age" )
  {
   Axis( -3000, -400);
   Boundary( "Start of Acton & Taunton");
   Phase( )
   {
    R_Date( "OxA-5948", 3225, 65);
    R_Date( "OxA-4651", 3220, 80);
    R_Date( "OxA-5949", 3110, 50);
    R_Date( "OxA-6177", 3055, 50);
    R_Date( "OxA-5196", 3035, 40);
    Interval( "Duration of Acton & Taunton");
   };
   Boundary( "Acton & Taunton - Penard");
   Phase( )
   {
    R_Date( "OxA-5187", 3045, 55);
    R_Date( "OxA-1526", 3030, 100);
    R_Date( "OxA-5953", 3015, 45);
    R_Date( "OxA-5951", 2980, 45);
    R_Date( "OxA-4504", 2965, 65);
    R_Date( "OxA-5952", 2965, 45);
    R_Date( "OxA-5959", 2965, 45);
    R_Date( "OxA-5183", 2930, 40);
    R_Date( "OxA-4653", 2910, 55);
    R_Date( "OxA-5950", 2910, 45);
    R_Date( "OxA-5954", 3025, 55);
    R_Date( "HAR-2940", 3020, 70);
    Interval( "Duration of Penard");
   };
   Boundary( "Penard - Wilburton");
   Page( );
   Axis( -3000, -400);
   Phase( )
   {
    R_Date( "OxA-4656", 3005, 75);
    R_Date( "OxA-4503", 2930, 55);
    R_Date( "OxA-4502", 2925, 50);
    R_Date( "OxA-5036", 2920, 50);
    R_Date( "OxA-5197", 2910, 50);
    R_Date( "OxA-5955", 2900, 45);
    R_Date( "OxA-5035", 2900, 45);
    R_Date( "OxA-5034", 2890, 45);
    R_Date( "OxA-5956", 2850, 50);
    R_Date( "OxA-5198", 2820, 70);
    Interval( "Duration of Wilburton");
   };
   Boundary( "Wilburton - Intermediate");
   Phase( )
   {
    R_Date( "OxA-5186", 2840, 40);
    R_Date( "OxA-5184", 2830, 65);
    R_Date( "OxA-5185", 2770, 50);
    Interval( "Duration of Intermediate");
   };
   Boundary( "Intermediate - Ewart Park");
   Page( );
   Axis( -3000, -400);
   Phase( )
   {
    R_Date( "OxA-5957", 2810, 45);
    R_Date( "OxA-4716", 2780, 50);
    R_Date( "OxA-4654", 2765, 45);
    R_Date( "OxA-5976", 2740, 45);
    R_Date( "OxA-4652", 2720, 45);
    R_Date( "BM-798", 2704, 50);
    R_Date( "OxA-5962", 2685, 60);
    R_Date( "OxA-6176", 2655, 50);
    R_Date( "OxA-5977", 2620, 45);
    Interval( "Duration of Ewart Park");
   };
   Boundary( "End of Ewart Park");
   Axis( -3000, -400);
  };

Note: these models will take some time to run so only do so if you are particularly interested in the results. You may, however wish to see how the models look under the different views of the project manager.

See also [Tools | Constraints | Groupings]


Estimating the end of a phase

Here we look at a couple of examples of estimating the end date of a phase. In this case we have a series of antler picks in a ring ditch (data is from Cleal et al. 1995). In the first model we assume that there was a period over which antlers were being collected and used and that the samples we have are likely to be a representative set of samples from this period.

  Sequence()
  {
   Boundary( "Start of antler collection");
   Phase( "Ditch antlers")
   {
    R_Date( "UB-3788", 4381, 18);
    R_Date( "UB-3787", 4375, 19);
    R_Date( "UB-3789", 4330, 18);
    R_Date( "UB-3790", 4367, 18);
    R_Date( "UB-3792", 4365, 18);
    R_Date( "UB-3793", 4393, 18);
    R_Date( "UB-3794", 4432, 22);
    R_Date( "BM-1583", 4410, 60);
    R_Date( "BM-1617", 4390, 60);
    Interval( "Period of antler collection");
   };
   Boundary( "End of antler collection");
  };

The date for the end of antler collection might be taken as the most likely date for the end of construction of the ditch. However, it is actually most likely that the antlers found come from the last stage of construction - with a few older antlers being mixed in the same context. This might be better modelled by an exponential distribution - rising to greatest concentration of samples found from the end of construction as in this model:

  Sequence()
  {
   Tau_Boundary( "T");
   Phase( "ditch antlers")
   {
    R_Date( "UB-3788", 4381, 18);
    R_Date( "UB-3787", 4375, 19);
    R_Date( "UB-3789", 4330, 18);
    R_Date( "UB-3790", 4367, 18);
    R_Date( "UB-3792", 4365, 18);
    R_Date( "UB-3793", 4393, 18);
    R_Date( "UB-3794", 4432, 22);
    R_Date( "BM-1583", 4410, 60);
    R_Date( "BM-1617", 4390, 60);
    Interval( "period ditch antlers");
   };
   Boundary( "E");
  };
  Tau=E-T;
  Tau&=U(0,50);

As an additional elaboration of the model in this case we have both calculated the time constant associated with the exponential (with the command Tau=E-T;) and we have also put a constraint on this, allowing it to be anything up to 60 years (with the command Tau&=U(0,50);). You would only do this in a model if you had reason for such an assumption.

See also [Tools | Groupings]


Wiggle matching tree ring sequences

Radiocarbon dating wood with visible tree-rings allows us to make good use of the wiggles in the calibration curve (see Bronk Ramsey et al. 2001). The easiest way to set up such models is to use the [Tools > Models > Dendro wiggle match] tool.

Typically you measure wood in decadal blocks. If the outer block contains the bark edge the felling date will be about five years later than the middle of the final block and so a model like this (data from Bronk Ramsey et al. 2004 or Manning et al. 2006) is suitable:

Options()
{
 Resolution=1;
};
Plot()
{
 D_Sequence( "Wiggle-match example")
 {
  R_Date( "P-14095", 3413, 22);
  Gap( 10);
  R_Date( "P-14096", 3430, 23);
  Gap( 10);
  R_Date( "P-14097", 3432, 22);
  Gap( 10);
  R_Date( "P-14098", 3431, 22);
  Gap( 10);
  R_Date( "P-14099", 3379, 22);
  Gap( 10);
  R_Date( "P-14100", 3371, 23);
  Gap( 10);
  R_Date( "P-14101", 3371, 22);
  Gap( 5);
  Date("Felling date");
 };
};

See also [Tools | Deposition models]


Sapwood estimates for dendro-chronology

OxCal provides methods for estimating the number of sapwood rings in cases where the heartwood-sapwood boundary is present but you do not have the bark edge. In order to use this facility you need to use a model based on actual sapwood measurements. Such a model has been worked out from data on mainland Britain (Miles 2005). This model is specifically recommended for post-Roman oak from England an Wales. Tools are provided for analysis of data from other regions (see Sapwood estimates).

In the following example you can see how results from several different timbers can be used together to provide a more precise date for the ceiling structure of the Jerusalem Chamber undercroft at Westminster Abbey. Note the setting of the resolution to 1 (the default of 5 is too coarse), the line specifying the model and then the use of the Combine() function to combine the results for the different samples.

 Options()
 {
  Resolution=1;
 };
 Plot()
 {
  Sapwood_Model("Mainland Britain", 2.77292, 0.100001, -0.275445, 0.314286377);
  Combine()
  {
   Sapwood("wa21", 1329, 243, 0, 1.06);
   Sapwood("wa22", 1354, 58, 6, 2.74);
   Sapwood("wa23", 1342, 55, 0, 2.55);
   Sapwood("wa26", 1328, 62, 0, 1.71);
   Sapwood("wa28a", 1353, 86, 0, 1.48);
   Sapwood("wa24a", 1337, 76, 0, 1.61);
  };
 };

Note that in this case the agreement is not quite as good as you would expect, possibly because of stockpiling of wood - so all of the trees involved might not have been felled in the same year. If you have an estimate for stockpiling times, these can be added into the calculation by adding the time between stockpiling and construction of the structure as in:

  Sapwood_Model("Mainland Britain", 2.77292, 0.100001, -0.275445,0.314286377);
  Combine()
  {
   Sapwood("wa21", 1329, 243, 0, 1.06)
    +Prior("Stockpile");
   Sapwood("wa22", 1354, 58, 6, 2.74)
    +Prior("Stockpile");
   Sapwood("wa23", 1342, 55, 0, 2.55)
    +Prior("Stockpile");
   Sapwood("wa26", 1328, 62, 0, 1.71)
    +Prior("Stockpile");
   Sapwood("wa28a", 1353, 86, 0, 1.48)
    +Prior("Stockpile");
   Sapwood("wa24a", 1337, 76, 0, 1.61)
    +Prior("Stockpile");
  };

In this case the combined date is for building the structure rather than for the felling of the trees. The file Stockpile.prior will need to contain probabilities of stockpiling for various numbers of years (see file formats).

See also [Sapwood estimates]


Deposition models

Deposition models can most easily be set up by using the [Tools > Models > Deposition models] dialogue. Here we show some of the deposition models that can be implemented with the program. In all cases we will assume the same radiocarbon dates for samples from a particular depth in a sequence. In the first model we will include the depth information but not use it directly. Here we use the simple Sequence() to specify that the samples were deposited in the sediment in the order that they are found. The depths are specified in meters (any unit can be used) and are added as attributes for each sample.

  Sequence("Simple sequence")
  {
   Boundary("Bottom"){};
   R_Date("",1010,25){ z=0.65; };
   R_Date("",887,25) { z=0.61; };
   R_Date("",979,25) { z=0.57; };
   R_Date("",848,25) { z=0.53; };
   R_Date("",809,25) { z=0.49; };
   R_Date("",743,25) { z=0.44; };
   R_Date("",595,25) { z=0.38; };
   R_Date("",613,25) { z=0.32; };
   R_Date("",485,25) { z=0.26; };
   R_Date("",395,25) { z=0.20; };
   Boundary("Top")   {};
  };

After running such a model we can plot the results against depth - but the depth information is not used in the calculations.

In the next model, we will make the assumption that the deposition rate is constant from 0.65m to 0.47m, it then changes and is constant again from 0.47m to 0.2m depth. This model is specified in terms of a U_Sequence().

  U_Sequence("Linear deposition")
  {
   Boundary("Bottom"){};
   R_Date("",1010,25){ z=0.65; };
   R_Date("",887,25) { z=0.61; };
   R_Date("",979,25) { z=0.57; };
   R_Date("",848,25) { z=0.53; };
   R_Date("",809,25) { z=0.49; };
   Boundary("Change"){ z=0.47; };
   R_Date("",743,25) { z=0.44; };
   R_Date("",595,25) { z=0.38; };
   R_Date("",613,25) { z=0.32; };
   R_Date("",485,25) { z=0.26; };
   R_Date("",395,25) { z=0.20; };
   Boundary("Top")   {};
  };

This model allows for no fluctuations in deposition rate. By using the P_Sequence() model instead we can allow for fluctuations. In the following example we assume that the deposition is characterised by an event size of 1 mm (that is the k value is 1000 events per meter):

  P_Sequence("Random variation",1000)
  {
   Boundary("Bottom"){};
   R_Date("",1010,25){ z=0.65; };
   R_Date("",887,25) { z=0.61; };
   R_Date("",979,25) { z=0.57; };
   R_Date("",848,25) { z=0.53; };
   R_Date("",809,25) { z=0.49; };
   Boundary("Change"){ z=0.47; };
   R_Date("",743,25) { z=0.44; };
   R_Date("",595,25) { z=0.38; };
   R_Date("",613,25) { z=0.32; };
   R_Date("",485,25) { z=0.26; };
   R_Date("",395,25) { z=0.20; };
   Boundary("Top")   {};
  };

Here the uncertainties are higher. If you reduce the k value to 500, the uncertainty will be higher still. Fine sediments would be expected to have higher k values than coarse ones.

Finally we assume that the underlying rate of deposition is not approximately constant over the two sections of the sequence but that it is exponentially increasing gradually over the whole range. We still allow for random fluctuation in the deposition rate, by using the P_Sequence() model.

  P_Sequence("Exponential rise",1000)
  {
   Tau_Boundary(""){};
   R_Date("",1010,25){ z=0.65; };
   R_Date("",887,25) { z=0.61; };
   R_Date("",979,25) { z=0.57; };
   R_Date("",848,25) { z=0.53; };
   R_Date("",809,25) { z=0.49; };
   R_Date("",743,25) { z=0.44; };
   R_Date("",595,25) { z=0.38; };
   R_Date("",613,25) { z=0.32; };
   R_Date("",485,25) { z=0.26; };
   R_Date("",395,25) { z=0.20; };
   Boundary("Top")   {};
  };

See also [Tools | Deposition models]