diff --git a/Manifest.toml b/Manifest.toml
index 4701490c43bf37648fa48a3db6c3b9024fd36546..07ccd8bde84b9b9208b8ddc15cdcd97407872a09 100644
--- a/Manifest.toml
+++ b/Manifest.toml
@@ -2,7 +2,7 @@
 
 julia_version = "1.9.3"
 manifest_format = "2.0"
-project_hash = "88b08cc01ff4cf4b3ac05aaa043f66221dec37b4"
+project_hash = "a93906f69837f39319eb070be15a9e813d23a9ed"
 
 [[deps.AbstractFFTs]]
 deps = ["ChainRulesCore", "LinearAlgebra"]
@@ -161,6 +161,12 @@ git-tree-sha1 = "ded953804d019afa9a3f98981d99b33e3db7b6da"
 uuid = "944b1d66-785c-5afd-91f1-9de20f533193"
 version = "0.7.0"
 
+[[deps.ColorSchemes]]
+deps = ["ColorTypes", "ColorVectorSpace", "Colors", "FixedPointNumbers", "PrecompileTools", "Random"]
+git-tree-sha1 = "67c1f244b991cad9b0aa4b7540fb758c2488b129"
+uuid = "35d6a980-a343-548e-a6ea-1d62b119f2f4"
+version = "3.24.0"
+
 [[deps.ColorTypes]]
 deps = ["FixedPointNumbers", "Random"]
 git-tree-sha1 = "eb7f0f8307f71fac7c606984ea5fb2817275d6e4"
@@ -347,9 +353,9 @@ version = "0.1.1"
 
 [[deps.FileIO]]
 deps = ["Pkg", "Requires", "UUIDs"]
-git-tree-sha1 = "7be5f99f7d15578798f338f5433b6c432ea8037b"
+git-tree-sha1 = "299dc33549f68299137e51e6d49a13b5b1da9673"
 uuid = "5789e2e9-d7fb-5bc7-8068-2c6fae9b9549"
-version = "1.16.0"
+version = "1.16.1"
 
 [[deps.FilePathsBase]]
 deps = ["Compat", "Dates", "Mmap", "Printf", "Test", "UUIDs"]
diff --git a/Project.toml b/Project.toml
index 91651ed6adbbe90e633fdad8d8b05045aa3f10a2..ed045a943924fb3cdee9caa015db6a7a165ddafc 100644
--- a/Project.toml
+++ b/Project.toml
@@ -7,15 +7,18 @@ version = "0.2.0"
 Agents = "46ada45e-f475-11e8-01d0-f70cc89e6671"
 ArgParse = "c7e460c6-2fb9-53a9-8c5b-16f535851c63"
 CSV = "336ed68f-0bac-5ca0-87d4-7b16caf5d00b"
+ColorSchemes = "35d6a980-a343-548e-a6ea-1d62b119f2f4"
 DataFrames = "a93c6f00-e57d-5684-b7b6-d8193f3e46c0"
 DataFramesMeta = "1313f7d8-7da2-5740-9ea0-a2ca25f37964"
 Dates = "ade2ca70-3891-5945-98fb-dc099432e06a"
 Distributed = "8ba89e20-285c-5b6f-9357-94700520ee1b"
+FileIO = "5789e2e9-d7fb-5bc7-8068-2c6fae9b9549"
 GeoArrays = "2fb1d81b-e6a0-5fc5-82e6-8e06903437ab"
 Logging = "56ddb016-857b-54e1-b83d-db4d58db5568"
 LoggingExtras = "e6f89c97-d47a-5376-807f-9c37f3926c36"
 Pkg = "44cfe95a-1eb2-52ea-b672-e2afdf69b78f"
 Random = "9a3f8284-a2c9-5f02-9a11-845980a1fd5c"
+Serialization = "9e88b42a-f829-5b0c-bbe9-9e923198166b"
 StableRNGs = "860ef19b-820b-49d6-a774-d7a799459cd3"
 StatsBase = "2913bbd2-ae8a-5f71-8c99-4fb6c76f3a91"
 TOML = "fa267f1f-6049-4f14-aa54-33bafae1ed76"
diff --git a/src/Persefone.jl b/src/Persefone.jl
index dd56cdba309bb88de6c5929c3bfb6132ad11d83f..165f3992f082c4aecbd85917456325d981e4cd1e 100644
--- a/src/Persefone.jl
+++ b/src/Persefone.jl
@@ -21,10 +21,13 @@ using
     DataFrames,
     DataFramesMeta,
     Distributed,
+    FileIO,
     GeoArrays, #XXX this is a big dependency - can we get rid of it?
     Logging,
     LoggingExtras,
+    #CairoMakie,
     Random,
+    Serialization,
     StableRNGs,
     TOML
 
@@ -76,12 +79,18 @@ export
     initialise,
     stepsimulation!,
     createevent!,
-    finalise!
+    finalise!,
+    #visualisemap,
+    #populationtrends,
+    #visualiseoutput,
+    savemodelobject,
+    loadmodelobject
 
 ## include all module files (note that the order matters - if file
 ## b references something from file a, it must be included later)
 include("core/input.jl")
 include("core/output.jl")
+#include("analysis/makieplots.jl")
 
 include("world/landscape.jl")
 include("world/weather.jl")
diff --git a/src/core/input.jl b/src/core/input.jl
index 7869c0e21e9bc5b1b32086b4606f612a5295dcf1..7d3ccf7d50c42e56f555701775754d70d80350ac 100644
--- a/src/core/input.jl
+++ b/src/core/input.jl
@@ -163,7 +163,7 @@ function parsecommandline()
             arg_type = String
             required = false
     end
-    #XXX this changes the global RNG?! (https://github.com/carlobaldassi/ArgParse.jl/issues/121)
+    #XXX this changes the global RNG?! (https://github.com/carlobaldassi/ArgParse.jl/issues/121) -> should be fixed in Julia 1.10
     args = parse_args(s)
     for a in keys(args)
         (args[a] == nothing) && delete!(args, a)
@@ -171,3 +171,24 @@ function parsecommandline()
     args
 end
 
+"""
+    loadmodelobject(fullfilename)
+
+Serialise a model object and save it to file for later reference.
+Includes the current model and Julia versions for compatibility checking.
+"""
+function loadmodelobject(fullfilename::String)
+    object = deserialize(fullfilename)
+
+    if !(typeof(object) <: Dict && typeof(object["model"]) <: AgentBasedModel)
+        @warn "This file does not contain a model object. Loading failed."
+        return
+    end
+    if (object["modelversion"] != pkgversion(Persefone) ||
+        object["juliaversion"] != VERSION)
+        @warn "This model object was saved with a different version of Persefone or Julia. It may be incompatible."
+    end
+    model = object["model"]
+    model.logger = modellogger(@param(core.loglevel), @param(core.outdir)) #reset logger
+    model
+end
diff --git a/src/core/output.jl b/src/core/output.jl
index debff8c7dfc93e4b78b30120eff07034328d8253..0ce93ea976709348b8f1552376367af6ae05c35e 100644
--- a/src/core/output.jl
+++ b/src/core/output.jl
@@ -20,7 +20,7 @@ function createdatadir(outdir::String, overwrite::Union{Bool,String})
             println("Type 'yes' to overwrite this directory. Otherwise, the simulation will abort.")
             print("Overwrite? ")
             answer = readline()
-            (answer == "yes") && (overwrite = true)
+            overwrite = (answer == "yes" || answer == "y")
         end
         !overwrite ? Base.error("Output directory exists, will not overwrite. Aborting.") :
             @warn "Overwriting existing output directory $(outdir)."
@@ -131,12 +131,14 @@ Struct fields:
     - header: a list of column names
     - outputfunction: a function that takes a model object and returns data values to record (formatted as a vector of vectors)
     - frequency: how often to call the output function (daily/monthly/yearly/end/never)
+    - plotfunction: a function that takes a model object and returns a Makie figure object (optional)
 """
 struct DataOutput
     name::String
     header::Vector{String}
     outputfunction::Function
     frequency::String
+    plotfunction::Union{Function,Nothing}
 end
 
 """
@@ -146,11 +148,12 @@ Create and register a new data output. This function must be called by all submo
 that want to have their output functions called regularly.
 """
 function newdataoutput!(model::AgentBasedModel, name::String, header::Vector{String},
-                        outputfunction::Function, frequency::String)
+                        outputfunction::Function, frequency::String,
+                        plotfunction::Union{Function,Nothing}=nothing)
     if !(frequency in ("daily", "monthly", "yearly", "end", "never"))
         Base.error("Invalid frequency '$frequency' for $name.") #TODO replace with exception
     end
-    ndo = DataOutput(name, header, outputfunction, frequency)
+    ndo = DataOutput(name, header, outputfunction, frequency, plotfunction)
     append!(model.dataoutputs, [ndo])
     if frequency != "never"
         if @param(core.csvoutput)
@@ -204,3 +207,36 @@ function outputdata(model::AgentBasedModel)
         end
     end
 end
+
+# """
+#     visualiseoutput(model)
+
+# Cycle through all data outputs and call their respective plot functions,
+# saving each figure as a PDF.
+# """
+# function visualiseoutput(model::AgentBasedModel)
+#     #TODO write tests
+#     @debug "Visualising output."
+#     for output in model.dataoutputs
+#         isnothing(output.plotfunction) && continue
+#         figure = output.plotfunction(model)
+#         save(joinpath(@param(core.outdir), output.name*".pdf"), figure)
+#     end
+# end
+
+"""
+    savemodelobject(model, filename)
+
+Serialise a model object and save it to file for later reference.
+Includes the current model and Julia versions for compatibility checking.
+
+WARNING: produces large files (>100 MB) and takes a while to execute.
+"""
+function savemodelobject(model::AgentBasedModel, filename::String)
+    object = Dict("model"=>model,
+                  "modelversion"=>pkgversion(Persefone),
+                  "juliaversion"=>VERSION)
+    filename = joinpath(@param(core.outdir), filename*".dat")
+    serialize(filename, object)
+    @debug "Saved model object to $(filename)."
+end
diff --git a/src/core/simulation.jl b/src/core/simulation.jl
index ea8a26e73a0dea6c370a21628ad69c6b633d8c6f..21bff1a5dcb0ce5fa996d15b3394f48b118d4c1a 100644
--- a/src/core/simulation.jl
+++ b/src/core/simulation.jl
@@ -26,6 +26,7 @@ end
 Carry out a complete simulation run using a pre-initialised model object.
 """
 function simulate!(model::AgentBasedModel)
+    @info "Simulation run started at $(Dates.now())."
     runtime = Dates.value(@param(core.enddate)-@param(core.startdate))+1
     step!(model, dummystep, stepsimulation!, runtime)
     finalise!(model)
@@ -42,7 +43,6 @@ parameters, in which case it creates a full-factorial simulation experiment
 and returns a vector of model objects.
 """
 function initialise(config::String=PARAMFILE, seed::Union{Int64,Nothing}=nothing)
-    @info "Simulation run started at $(Dates.now())."
     settings = getsettings(config, seed)
     scanparams = settings["internal.scanparams"]
     delete!(settings, "internal.scanparams")
@@ -152,7 +152,7 @@ function finalise!(model::AgentBasedModel)
     with_logger(model.logger) do
         @info "Simulated $(model.date-@param(core.startdate))."
         @info "Simulation run completed at $(Dates.now()),\nwrote output to $(@param(core.outdir))."
-        #XXX is there anything to do here?
+        #visualiseoutput(model) #TODO
         model
     end
 end
diff --git a/src/world/landscape.jl b/src/world/landscape.jl
index 1cacf092762467096bddead5430d539c42aa52ad..d97a91c01ee3a3da849acaa3791da42f02e40aa4 100644
--- a/src/world/landscape.jl
+++ b/src/world/landscape.jl
@@ -43,9 +43,11 @@ configuration. Returns a matrix of pixels.
 """
 function initlandscape(landcovermap::String, farmfieldsmap::String)
     @debug "Initialising landscape"
+    #TODO replace errors with exception
+    !(isfile(landcovermap)) && Base.error("Landcover map $(landcovermap) doesn't exist.")
+    !(isfile(farmfieldsmap)) && Base.error("Field map $(farmfieldsmap) doesn't exist.")
     landcover = GeoArrays.read(landcovermap)
     farmfields = GeoArrays.read(farmfieldsmap)
-    #TODO replace error with exception
     (size(landcover) != size(farmfields)) && Base.error("Input map sizes don't match.")
     width, height = size(landcover)[1:2]
     landscape = Matrix{Pixel}(undef, width, height)
diff --git a/test/io_tests.jl b/test/io_tests.jl
index 435bd253c25385511cce3c344484d4fafde8cdc0..7605365e3243dfe4f1a53b6d96bfbedd27743be3 100644
--- a/test/io_tests.jl
+++ b/test/io_tests.jl
@@ -63,3 +63,21 @@ end
     @test size(model.datatables["end"]) == (1, 1)
     rm(outdir, force=true, recursive=true)
 end
+
+@testset "Model object serialization" begin
+    model = inittestmodel()
+    Ps.createdatadir(@param(core.outdir), @param(core.overwrite))
+    @test_logs((:debug, "Saved model object to results_testsuite/test.dat."),
+               min_level=Logging.Debug, match_mode=:any,
+               savemodelobject(model, "test"))
+    @test isfile(joinpath(@param(core.outdir), "test.dat"))
+    model2 = loadmodelobject(joinpath(@param(core.outdir), "test.dat"))
+    @test model.date == model2.date
+    @test model.settings == model2.settings
+    @test length(model.agents) == length(model2.agents)
+    simulate!(model)
+    simulate!(model2)
+    @test model.date == model2.date
+    @test length(model.agents) == length(model2.agents)
+    rm(@param(core.outdir), force=true, recursive=true)
+end
diff --git a/test/runtests.jl b/test/runtests.jl
index dc86b72f29752292c7427b6868526835af8ede2a..f2dca73a018944e26ae461def5c1558ae988476e 100644
--- a/test/runtests.jl
+++ b/test/runtests.jl
@@ -22,6 +22,8 @@ const TESTPARAMETERS = joinpath(pkgdir(Persefone), "test/test_parameters.toml")
 const TESTSETTINGS = Ps.getsettings(TESTPARAMETERS)
 
 """
+    inittestmodel(smallmap=true)
+
 Initialise an AgentBasedModel for testing purposes.
 
 `smallmap`: use a hypothetical small landscape rather than a real one?