Widespread unexplained variations in clinical practices and patient outcomes, together with rapidly growing availability of data, suggest major opportunities for improving the quality of medical care. One way that healthcare practitioners try to do that is by participating in organized healthcare quality improvement collaboratives (QICs). In QICs, teams of practitioners from different hospitals exchange information on clinical practices, with the aim of improving health outcomes at their own institutions. However, what works in one hospital may not work in others with different local contexts, due to non-linear interactions among various demographics, treatments, and practices. I.e., the clinical landscape is a complex socio-technical system that is difficult to search. In this dissertation we develop methods for analysis and modeling of complex systems, and apply them to the problem of healthcare improvement. Searching clinical landscapes is a multi-objective dynamic problem, as hospitals simultaneously optimize for multiple patient outcomes. We first discuss a general method we developed for finding which changes in features may be associated with various changes in outcomes at different points in time with different delays in affect. This method correctly inferred interactions on synthetic data, however the complexity and incompleteness of the real hospital dataset available to us limited the usefulness of this approach. We then discuss an agent-based model (ABM) of QICs to show that teams comprising individuals from similar institutions outperform those from more diverse institutions, under nearly all conditions, and that this advantage increases with the complexity of the landscape and the level of noise in assessing performance. We present data from a network of real hospitals that provides encouraging evidence of a high degree of similarity in clinical practices among hospitals working together in QIC teams. Based on model outcomes, we propose a secure virtual collaboration system that would allow hospitals to efficiently identify potentially better practices in use at other institutions similar to theirs, without any institutions having to sacrifice the privacy of their own data. To model the search for quality improvement in clinical fitness landscapes, we need benchmark landscapes with tunable feature interactions. NK landscapes have been the classic benchmarks for modeling landscapes with epistatic interactions, but the ruggedness is only tunable in discrete jumps. Walsh polynomials are more finely tunable than NK landscapes, but are only defined on binary alphabets and, in general, have unknown global maximum and minimum. We define a different subset of interaction models that we dub as NM landscapes. NM landscapes are shown to have smoothly tunable ruggedness and difficulty and known location and value of global maxima. With additional constraints, we can also determine the location and value of the global minima. The proposed NM landscapes can be used with alphabets of any arity, from binary to real-valued, without changing the complexity of the landscape. NM landscapes are thus useful models for simulating clinical landscapes with binary or real decision variables and varying number of interactions. NM landscapes permit proper normalization of fitnesses so that search results can be fairly averaged over different random landscapes with the same parameters, and fairly compared between landscapes with different parameters. In future work we plan to use NM landscapes as benchmarks for testing various algorithms that can discover epistatic interactions in real world datasets.